Methods and systems for drug screening and computational modeling based on biologically realistic neurons

ABSTRACT

A method for screening a test composition for potential efficacy in treatment of a disorder includes a first computer model representative of a volume of disease-afflicted neural tissue comprising biologically realistic neurons exposed to the test composition; and providing an initial excitation to the first computer model. Following a selected computation interval, a first outcome is determined. The first outcome indicates a response of the first computer model to the initial excitation and indicates whether the test composition has the potential to be effective in treating the disorder.

RELATED APPLICATIONS

Under 35 USC 119, this application claims the benefit of the prioritydate of U.S. provisional application 60/736,599, filed on Nov. 14, 2005,the contents of which are herein incorporated by reference.

FIELD OF THE INVENTION

This invention relates to methods and systems for screening potentialmedications and compositions using biologically realistic computationalnetwork models. The computational network model can be used as a moderndrug discovery tool for a variety of medical disorders.

BACKGROUND OF THE INVENTION

The human brain is made up of neurons connected to one another in acomplex network. It is believed that when humans learn, new connectionsare made or existing ones are modified. Neural networks, used inartificial intelligence (A.I.) applications, are massively parallelcomputing models inspired by the human brain. Such networks aretypically implemented by multiple processors connected by adaptiveweights. Computational neural models based on such neural networks cansimulate brain conditions and provide valuable information and insightabout the human brain.

Doctors who treat patients having neuropsychiatric disorders often relyon psychiatric drugs. Specific treatments and drugs exist for specificdiagnostic categories of patients. For example, neuroleptics areprescribed for schizophrenia, antidepressants are administered fordepression, anxiolytics for anxiety, lithium for mania, and stimulants,such as RITALIN®, for attention-deficit/hyperactivity disorder (ADHA).

Before prescribing drugs to humans, it is prudent to first test them toensure that they are effective, or at the very least, that they aresafe. The most straightforward way to test such drugs is to test them onhumans. As late as the middle of the last century, for example, it wasroutine to test drugs on prisoners and on patients in mental asylums. Infact, the well-known antipsychotic drug chlorpromazine was discovered inthis way.

Since then, testing new drugs in this way has been acknowledged asunethical. As a result, the safety and efficacy of new drugs is assessedby testing them on animal models. In the case of psychiatric drugs, thisgenerally involves isolating an animal behavior that is thought to beanalogous to a human behavior or psychiatric condition, administeringthe drug in question to the animal, and observing if the behaviorchanges. As an example, the operant conflict test in rate is thought toembody behavior analogous to anxiety in humans and the medicationDIAZEPAM®, also known as VALIUM®, was discovered by its ability todecrease such behavior in rats in a laboratory environment.

However, not all psychiatric disorders and neuropsychiatric conditionshave clear behavioral correlates in animals. For certainneuropsychiatric brain disorders, such as schizophrenia, it is difficultto test the efficacy of a potential drug preclinically, in part becausethere are presently no well-established animal models of schizophrenia.As a result, using animal models to screen antipsychotic medications forantischizophrenic potency is not feasible. Thus, for many neurologicandisorders, there are currently no well-understood and establishedmethods for preclinically screening potential drugs safely andexpeditiously.

SUMMARY

The invention features screening test compistitions for their efficacyin the treatment of various types of medical, e.g., neural disordersprior to human clinical trials. Such neural disorders can include, forexample, brain disorders manifested as psychiatric disorders, or otherdisorders, for example disorders of motor function, that are susceptibleto neurocomputational modeling.

In one aspect, the invention includes a method for screening a testcomposition for potential efficacy in treatment of a disorder (e.g.,schizophernia). The method includes creating a first computer model(e.g., an “attractor” neural network model) representative of a volumeof disease-afflicted neural tissue (e.g., hippocampus) exposed to thetest composition. The first computer model includes interconnectedneurons, in each of which the quantitatively simulates spatialiycompartmentalized ion channel conductance and receptor actvity. Aninitiai excitation (e.g., a biologically realistic one) is then providedto the first computer model, and an outcome is then determined. Thisfirst outcome is indicative of a response of the first computer model tothe initial excitation, which also represents the potential efficacy ofthe test composition for treating the disorder.

The method can include modeling neurons having an anatomically realisticdendritic arborization (e.g., based on anatomical data from actualcells). The method can also include modeling neurons that have varyingion channel and receptor distributions along the dendro-somatic axis(e.g., based on experimental data for ion channel or receptordistribution). The first computer model can also simulate changes insecond messenger concentrations according to the varying distributionsof ion channels and receptors. In addition, the first computer model canalso simulate the presence ion channel or receptor antagonists.Likewise, the model can also simulate the presence of ion channel orreceptor antagonists. In one embodiment, the model imports data from adatabase (e.g., data for: dendritic arborizations, ion channeldistribution, receptor distribution, ion channel antagonists, ionchannel agonists, receptor antagonists, or receptor agonists).

One practice of the invention includes creating a second computer modelrepresentative of the volume of disease-afflicted neural tissue. Thisvolume is not exposed to the test composition. An initial excitation isthen applied to this second computer model to produce a second outcomeindicative of a response of the second computer model to the initialexcitation. A difference between the first and second outcomes indicatesthat the test composition is a candidate composition for treating thedisorder.

In another practice of the invention, a third computer modelrepresenting a volume of neural tissue free of the disease is createdand provided with an initial excitation to produce a third outcomeindicative of a response of the third computer model to the initialexcitation. A similarity between the first and third outcomes indicatesthat the test composition is a candidate composition for treating thedisorder. In one embodiment, a comparison between the outcomes of afirst and third model includes determining a first epoch frequencyassociated with the first outcome, determining a test epoch frequencyassociated with the first outcome, outcomes, and classifying the testcomposition as a candidate composition for treating the disorder whenthe first epoch frequency is associated with the third outcome. (Thetest composition is also classified as a candidate for treating thedisorder when the first epoch frequency is greater than the test epochfrequency and the test epoch frequency is associated with the secondoutcome.

In yet another practice of the invention, creating a first computermodel of the disorder (e.g., a neuropsychiatric or neurologicaldisorder) includes defining a model of a volume of neural tissue (e.g.,defining a neural network), the model having a population profile ofbiologically realistic neurons. The model is then altered to simulatelesioning caused by the disorder (e.g., by altering the populationprofile to be consistent with the disorder) as well as to simulate theeffect of the test compisition on the neural tissue. In one embodiment,interneuron density is reduced in the first computer model, to alter thepopulation profile to be consistent with the disorder.

In another aspect, the invention includes a method for generating afirst computer model for determining whether a test composition is acandidate drug for treating a disorder (e.g., a neuropsychiatric orneurological disorder). Generating the first computer model includesgenerating a computational network model of a normal human brainmanifesting a plurality of normal characteristics of human behavior, andintroducing a digital representation of one or more physiologicallesions into the normal computational network model. The physiologicallesions are selected to be consistent with the suspected neuropathologyof the disorder. The method further includes applying pharmacologicaldata of the test composition to the computational network, and based onthe application of the pharmacological data determining the efficacy ofthe test composition for treating the disorder. A favorable outcome inthe network model indicates that the test composition is a candidatedrug to treat the disorder. As used herein, “favorable” means that thecomputer simulation, i.e., the computational model of the braindisorder, shows signs of physiological improvement consistent withimprovements in the brain disorder.

In another aspect of the invention, creating a computer model includescreating a first computational network model for the disorder (e.g., aneuropsychiatric or neurological disorder) and applying input data of atest composition to the first network model to obtain the resultingdata. The resulting data from the first network model are then comparedwith resulting data from a second network model simulating exposure to acomposition known to be effective for treating the disorder. Theefficacy of the test composition is then determined based on thecomparison between the resulting data from the first and second modelnetworks.

In one embodiment, generating the computational network model for thedisorder includes generating a model that manifests a plurality ofnormal characteristics of human behavior, and then introducing one ormore physiological lesions to the normal computational network modelconsistent with suspected neuropathy of the disorder. In anotherembodiment, applying input data of the test composition includessimulating dopamine induced effects on the computational network model.

Various aspects of the invention can include one or more of thefollowing features. For example, functional characteristics can be addedto the comnputational network model to generate a model of the medicaldisorder by degrading neurons of the computational network model in amanner analogous to that which occurs in humans afflicted with theneuropsychiatric or neurological disorder. The method can also includedetermining the efficacy of the test composition for treatment of thedisorder includes determining whether the application of the testcomposition modifies behaviors attributable to the disorder in abeneficial way. Applying the pharmacological data of the testcomposition can include modeling effects of the test composition ondendritic input integrating for producing an axonal output, or affectingneurotransmitter release properties of neurons in the computationalnetwork model. The medical disorder can be selected to be schizophrenia,Alzeheimer's disease, dementia, or a seizure disorder.

In another aspect, the invention includes a system for screening a testcomposition. The system includes a processor and a memory coupled to theprocessor. Encoded in the memory is software that, when executed, causesthe processor to generate a computational network model includingbiologically realistic neurons, and manifesting a neuropsychiatric orneurological disorder, to apply pharmacological data of the testcomposition to the computational network model; and to determine, basedon the application of the physiological data to the network model of theneuropsychiatric or neurological disorder, the efficacy of the testcomposition for treatment of the disorder. The biologically realisticneurons of the model can feature one or more of: anatomically realisticdendritic arborizations (e.g., based on anatomical data from actualcells), spatially compartmentalized ion channel conductances, spatiallycompartmentalized changes in second messenger concentrations, ionchannel or receptor distributions that vary along thedendrosomaticaxonal axis. In one embodiment, the system imports datainto the model from a database.

In one aspect of the invention, the software when executed, causes theprocessor to generate a computational network model manifesting a numberof normal characteristics of human behavior. The software also causesthe processor to introduce structural or functional lesions to thenormal computational network consistent with suspected neuropathology ofthe neuropsychiatric or neurological disorder. For example, the softwarecan cause the processor to add functional characteristics to generate amodel of the medical disorder by degrading neurons of the computationalnetwork model in a manner analogous to the dysfunction of neurons inhumans afflicted with the disorder. In yet another embodiment, thesoftware causes the processor to apply the pharmacological data of thetest composition on ion channels.

As an example of the above, the invention can include acomputer-readable medium on which is encoded a data structurerepresentative of a biologically-realistic model of a volume ofhippocampal tissue. The data structure includes data representative ofpopulation of biologically realistic neurons in each layer of thehippocampus, data representative of types of neurons in each layer ofthe hippocampus; and data representative of synaptic connections betweenneurons in the hippocampus.

In yet another aspect, the invention includes a method for determiningthe potential efficacy of a test composition for treating a neuoraldisorder. The method includes creating a first computer modelrepresentatice of a disease-afflicted neuron exposed to the testcomposition, in which the model quantitatively simulates spatiallycompartmentalized conductance of number of ion channels and receptors,as well as the effect of the test composition on the conductance of theion channels. The model is then provided with an initial excitation, anda first outcome is determined, which is indicatice of the response ofthe first computer model to the initial excitation. The first outcome isrepresentative of the ability of the composition to treat the neuraldisorder.

In one practice of the invention, the method includes creating a secondcomputer model representative of the disease-afflicted neuron isolatedfrom the test composition. The second computer model is then providedwith an initial exitation, and a second outcome is determined, which isindicatice of a response of the second computer model to the initialexcitation. A difference between the outcome of the first computer modeland the outcome of the second computer model to an initial excitation,indicates that the test composition is a candidate drug for treating thedisorder.

In another practice of the invention, the method includes creating athird computer model representatice of a neuron free of the disease. Thethird computer model is then provided with an initial excitation toproduce a third outcome indicative of a response of the third computermodel to the initial excitation. A similarity between the first andthird outcomes indicates that the test composition is a candidatecomposition for treating the disorder.

In various embodiments, the method for determining the potentialefficacy of a test composition for treating a neural disorder caninclude one or more of the following features. For example, creating afirst or second computer model of a neuron can include simulating: ionchannel or receptor distributions that vary along the dendro-somaticaxisbased on experimental data, anatomically realistic dendriticarborizations (e.g., based on data from actual cells), changes in secondmessenger concentrations, and the presence of an ion channel antagonist,ion channel agonist, receptor antagonist, or receptor agonist. In oneembodiment, the method includes importing data into the model from thedatabase.

The invention provides several advantages.

The computational network modeling based on neural networks for drugscreening provides rapid screening of chemical compositions in silicofor various psychiatric and neurological disorders, even where nosuitable animal models exist.

The present computational network modeling systems and methods can beused, for example, by pharmaceutical and biotechnology companies, aspart of a comprehensive drug discovery and development process, as wellas by academic and other laboratories studying different pharmacologicalaspects of neuropsychiatric drugs. The new systems and methods of thepresent invention can be used, alone or in combination with existingmethods, to improve the accuracy of the drug screening process in thescreening of drug compositions prior to clinical trials. As a result,drug compositions with a likelihood of clinical effectiveness can bereadily and quickly identified. This would increase the number of drugs,medications, and agents screened while at the same time decreasing theoverall costs of the screening process. Further, drugs and medicationsthat are more effective than currently available drugs can be madesafely available for treatment.

The computational network models also provide the ability to modelmultiple pathways or multiple molecular targets in a single assay.Additionally, the computational models provide the ability to assesstest compositions based on their effects on electrophysiology, not justbiochemistry, both at the singe cell and network level.

Unless otherwise defined, all technical and scientific terms used hereinhave the same meaning as commonly understood by one of ordinary skill inthe art to wich this invention belongs. Although methods and materialssimilar or equivalent to those described herein can be used in thepractice or testing of the present invention, suitable methods andmaterials are described below. All publications, patent applications,patents, and other references mentioned herein are incorporated byreference in their entirety. In case of conflict, the presentspecification, including definitions, will control. In addition, thematerials, methods, and examples are illustrative only and not intendedto be limiting.

Other features and advantages of the invention will be apparent from thefollowing detailed description, and from the claims.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a schematical diagram for developing and testing acomputational network model.

FIG. 2 is a schematic representation of a two-memory state for anine-unit atractor network.

FIG. 3 is a schematic representation of an energy landscape of anattractor network.

FIG. 4 is a schematic representation of a computational network of anattractor network.

FIG. 5 is a schematic of representation of the process for generating amodel of biologically realistic neuron.

FIG. 6 shows ion channel density data for four different ion channels(Na⁺, I_(KA), I_(h), Ca²⁺) for a number of different types of neuron.

FIG. 7 shows an example of a schematic representation of a cell and itscompartments in the NeuronDB database.

FIG. 8 show a schematic of a biologically realistic neuron created bythe program “makeproto.c” (top panel) which is based on a real neuron(bottom panel).

FIG. 9 shows a schematic representation of a network of intracellularinteractions which shows examples of receptor→second messenger→ionchannel effects in hippocampus CA1 neurons.

FIG. 10 shows a graphical representation of an example of secondmessenger biosynthetic enzyme activity as a function of receptor agonistconcentration.

FIG. 11 shows a graphical representation of the operational model ofdrug action.

FIG. 12 is a detailed flow diagram of the computational drug screeningor discovery process.

FIG. 13 is a flow diagram for applying potential drug effects to thecomputational model of FIG. 12.

FIG. 14 is a flow diagram for generating a computational network modelof a disorder in the computational drug screening process of FIG. 12.

FIG. 15 is a diagram model of the parameters that can be chaned forsimulating the effect of a drug.

FIG. 16 is a flow diagram of another computational drug screening ordiscovery process.

FIG. 17 is a computer system used to carry out the computational testscreening process of FIG. 12.

FIG. 18 is a software module used in the computer system of FIG. 17.

FIG. 19 is a graph indicating axonal and dendritic arborizations ofinterneuron species.

FIGS. 20A and 20B are graphs comparing modeled and measured spike ratesfor interneurons.

FIGS. 20C and 20D are graphs comparing modeled and measured spike ratesfor interneurons.

FIGS. 21A and 21B are graphs comparing modeled and measured inter-spikeintervals for pyramidal cells.

FIGS. 21C and 21D are graphs comparing modeled and measured inter-spikeintervals for interneurons.

FIG. 22A is a graph of the output of a virtual electrode in a model of anormal hippocampus.

FIG. 22B is a graph of a firing pattern of neurons from a model of anormal hippocampus.

FIG. 23 is a graph of a neuron firing pattern from a model of ahippocampus with I_(AHP) channels removed.

FIG. 24 is a graph of neuron firing patterns from a model ofschizophrenic hippocampus.

FIG. 25 is a graph of neuron firing patterns from the hippocampus modelof FIG. 24 after having been treated with a psychiatric drug.

FIG. 26 shows the first 25 lines (and header) from a file in the “.swc”format, which specifies dendritic morphology of an actual neuron.

FIG. 27 is a printout of the file “chan_dens,” which lists in tabularform the densities of various species of ion channels as a function ofdistance from the soma.

DETAILED DESCRIPTION

The new drug screening computational network model enables theidentification of compositions having desired pharmacological effects onmedical disorders. The methods and systems provide valuable infomationand research and development support for pharmaceutical and othercompanies by enabling sreening of test compositions in silico formedical disorders, including those lacking adequate animal models.Therefore, the new methods and systems provide an adjunct tohigh-throughput screening (HTS) methods and other tools of modern drugdiscovery for various types of medical disorders, e.g., neuropsychiatricdisorders such as schizophrenia.

First, the overall systems and methods of screening test compositons aredescribed. In later sections, the details of individual steps andelements of the new systems and methods are presented. Then,representative examples are set forth.

General Methodology

The new methods and systems include modeling, on a computer, the effectsthat different chemical compositions (e.g., drugs) or biological factoshave on specific neural disrders. Whan applied to disorders for which noappropriate animal model exists, the new methods and systems ofcomputer-based drug modeling can identify the molecules with a greaterlikelihood of being efficacious. The method can be used alone, or inconjunction with medical chemistry approaches, rational drug design, insilico computer-aided drug design, or HTS methods.

For example, in the case of disrders for whicfh a computational modelcan be constructed, medication effects can be applied at the level ofsingle cells. This is generally done by sifting the model parametersthat control the neurophysiologic behavior of the model's constituentneurons. Such model parameters can include those representing ionchannel conductance or synaptic response proterties. The effects ofthese changes on the overall behavior of the network can then beobserved. The computational network model thus translates the effects ofa drug on individual neurons into the effects of the drug on abiological system made up of many neurons.

In its most general sense, a computational network model includes anumber of interconnected neurons whose emergent, collective behavior isintended to emulate that of a biological system. The network cellularlevel parameters present in a network model and measured in a biologicalsystem can be of any type. Neural models range in complexity from simple“computational units” that sum a weighted vector of inputs to produceand output, to complex compartmental models that take into accounttime-varying conductances and Hodgkins-Huxley ion channels. In no caseare all underlying neurobiological variables and mechanisms known.Consequently, it is often necessary to begin with a hypothesis toexplain an unknown mechanism. This hypothesis is then refined as moreexperimental data becomes available.

The development and testing of a computational model for simulating abiological system is by nature an iterative process. Referring to FIG.1, that process begins with the formulation of a hypothesis (step 102)and the construction of a network model (step 104) embodying thathypothesis. The assumptions about the behavior of the biological systemsthat underlie the hypothesis are embodied in the equations and theparameters of the model and generally implemented on a computer.

Execution of the model (step 104) results in an output 106 that attemptsto simulate the behavior of a biological system. This output 106 is thenquantitatively compared with corresponding experimentally observedvalues (step 110) obtained from the actual biological system (step 108)to determine how well the output (step 106) matches the correspondingobserved values from the biological system. If the difference ormismatch is below some threshold, the model is considered suitable foruse in making predictions about the behavior of the biological system(step 112). Otherwise, the parameters of the model parameters (step104). This process continues until model outputs are within reasonableagreement with experimental values.

Neural Networks

A neural network is a computational network model composed of neuronsand connections between these neurons. The strength of each connectioncan be expressed as a “weight”. The activation of a given neuron isbased on the activation of the neurons that have connections directed atthat neuron and the weights associated with those connections.

Computer based neural network simulation, which was originally inspiredby artificial intelligence (A.I.) research, has been used to study humancognition, including psychiatric and neurological mental disorders.Computational network models to study human cognition have been modifiedto better emulate the dynamic nature of human cognitive processes.Accordingly, neural network models are formulated to capture theemergent properties of assemblies of neurons functioning together. Aswill be described further below, computational network models can beused to study normal and pathological human brain functions. Therefore,by using neural network simulation, human thought processes can bemolded by applying principals governing the dynamic interactions withinneural ensembles to gain a better understanding of their gross behaviorviewed as a cognitive process.

Attractor Networks

It is believed that thoughts, or cognitions, are represented in thebrain as patterns of neural activation. The attractor, or “Hopfield,”network formulation embodies the maner in which cognitions may beencoded and recalled at the level of ensembles of neurons. This networkformulation is useful, as it describes the manner in which clinicallyobservable phenomeno (e.g., memory tasks) can be operationalizedcomputationally.

An attractor network is generally conceptualized as an array of neurons,with each neuron connected to every other. Each neuron calculates theweighted sum of its inputs and applies a transfer function to determinewhether the neuron will be active or inactive. A “memory” consists of apattern of activation of the consituent neurons. One teaches the networkby activating the neurons of a particular memory and then applying aHebbian learning rule (i.e., interneuronal synaptic connections are notstatic, but can rather vary in strength, depending on the past activityof the constituent neurons), which adjusts the weighting of theinter-neuronal connections. A well-taught network, when presented with afragment of a memory, will return the complete, intact memory. Theexample presented below is implemented using an attractor network.

Attractor networks have three interrelated properties: (1) memories arestored in a content addressable manner, (2) they are represented in adistributed, rather than localized, form; and (3) the systems arecapable of generalization. These characteristics capture the essentialways that attractor networks mirror actual neurobiological processes.

In a content-addressable system, the memory does not “exist” in anyparticular location, but is spread over several processing units. Forexample, humans are able to recall items from memory based on a partialdescription of their contents, and can even do so if the description isnot entirely correct. In such a system, referred to as a “contentaddressable system,” a fragment of a memory provides access to thecomplete, stored memory. In contrast to the mannier in which computersassign each memory to a particular location, artificial neural networkshave the ability to recall complete patterns from fragmentary inputstimuli.

In a distributed representation, several neurons are involved in thestorage of a given memory, and a particular neuron will often be part ofthe brain's representation of several distinct memories. In part forthis reason, removal of one or a small number of neurons from a neuralnetwork does not excise a particular memory, nor does it cause asignificant decline in overall performance.

In a generalized system, if a network learns new information about anitem, that new information is also applied to similar memories stored inthe network. This is because memories are distributed among nodes of thenetwork. Hence, two similar memoreis may have a large number ofactivated nodes in common; that is, their patterns of activation willoverlap.

For example, referring to FIG. 2, in an attractor network 300, nineneurons are arranged in a 3×3 array, with each neuron reciprocallyconnected to every other. In a given cycle, each neuron adjusts itsinternal state based on the inputs from the other eight neurons in thenetwork 300. The activation state of each neuron will be designatedμ_(i). An active neuron is represented as +1 and an inactive neuron as−1. In FIG. 2, two memory states 302 and 304 to be stored are shown.

The synaptic strength of the axonal connection from neuron i to neuronj, designated as T_(i→j), measures the effect of neuron i on neuron j.This synaptic strength can be positive (excitatory), negative(inhibitory), or 0 (i.e., no axonal connection). The overall effect onneuron j from all other neurons in the system, designated S_(j) below,is given by:

${S_{j} = {\sum\limits_{i = 1}^{9}\;{T_{i - j} \times \mu_{i}}}},{i \neq {j.}}$

If this sum positive, the neuron is activated (i.e., set to +1); if itis negative the neuron is rendered inactive (i.e., set to −1). For agiven cycle, the activation state of each neuron in the system iscalculated in this way.

The first step in implementing an attractor network is to storememories. This is done by examining each pair of neurons in the system.For instance, if, for a given memory, the activation states of twoneurons are the same (i.e., +1 and −1 and −1), the synaptic strengthbetween them is increased. If they are different (+1 and −1), thesynaptic strength between them is decreased. This is carried out for allneurons, across all memories. Mathematically, with M memories, this canbe represented by a weight T:

${T_{i - j} = {\sum\limits_{m = 1}^{M}{\mu_{i}^{m} \times \mu_{j}^{m}}}},$

If two neurons in the system are both active (or both inactive) inseveral memories, the weight T connecting them is large. Conversely, ifthe activation levels of two particular neurons are different, theweight T is negative. Once the neural network has been trained, the setof weights associated with the network, i.e., the weight matrix,represents the learning of the network. With the formulation given here,the number of memories a network remembers is approximately 0.15 timesthe number of neurons in the system.

For memories 302 and 304 of FIG. 2, the inter-neuronal weights for thefirst few combinations is as follows:T ¹⁻¹=(−1×−1)+(1×1)=2T ²⁻¹=(−1×−1)+(1×1)=2T ³⁻¹=(1×−1)+(1×1)=0T ⁴⁻¹=(1×−1)+(−1×1)=−2

This calculation is carried out for all 9×9 neuronal connections. Thelearning rule results in symmetric weights: T_(1→j)=T_(j→1).

To test the performance of a network after learning, a fragment of oneof the stored memories is supplied to the network. For example, the topthree neurons of a start-up array can match those of memory 302:−1 −1 10 0 00 0 0

Next the weighted input to each neuron in the system is computed basedon this startup array. For neuron 1, S would be calculated as follows:

$\begin{matrix}{S_{i} = {\left( {T_{11} \times \mu_{1}} \right) + \left( {T_{21} \times \mu_{2}} \right) + \left( {T_{31} \times \mu_{3}} \right) + \ldots}} \\{= {\left( {2 \times {- 1}} \right) + \left( {2 \times {- 1}} \right) + \left( {0 \times 1} \right) + \ldots}} \\{= {- 4}}\end{matrix}$

Performing the above computation for each neuron in the system yieldsthe following set of S_(j)s:S₁ S₂ S₃ −4 −3 2S₄ S₅ S₆ =4 4 −2S₇ S₈ S₉ −4 2 −2

Applying threshold functions (if S_(ƒ)>0, neuron j is on, otherwise,neuron j is off), the following pattern emerges:−1 −1 11 1−1−1 1 −1

In attractor networks, an energy function is defined as follows:

$E = {{{- 1}/2}{\sum\limits_{i = 1}^{N}{\sum\limits_{j = 1}^{N}{T_{ij}\mu_{i}\mu_{j}}}}}$

Referring to FIG. 3, the energy function can be viewed as an energylandscape 400 having an energy level (z axis) and a network state-space(the x-y plane). The landscape 400 includes all 2^(N) possiblecombinations of the network (where N is the number of neurons in thenetwork). Each of the local minima of the network corresponds to amemory. As the network evolves, i.e., as the activation states of itsneurons are cyclically updated, the energy decreases. Thus, as thenetwork goes through successive cycles of updating, the state of thenetwork flows down the valleys of the energy function, analogous to aball rolling down hill. The network eventually stabilizes at one ofthese local minima, or “attractors.” Referring to FIG. 4, a schematic500 of the state-space of a model with three attractors represents a“projection” of the energy landscape onto the x-y plane. The state-spaceis divided into “basins of attraction.” In a neural network, if thememory cue is within one of these basins, it is drawn to the indicatedmemory state.

The change in the energy level of the network due to the change in stateof a given neuron μ_(i) (i.e., from −1 to +1, or from +1 to −1) is givenby:

${{\Delta\; E_{i}} = {{{- 1}/2}\;\Delta\;\mu_{i}{\sum\limits_{j = 1}^{N}{T_{ij}\mu_{j}}}}},\mspace{20mu}{j \neq i}$Biologically Realistic Computer Model of a Neuron

To prepare a computational model of a neural tissue sample, on firstcreates a biologically realistic computational model of a single neuron.As used herein, a biologically realistic model of a neuron is one thatexplicitly models trans-membrane potentials and time-variations thereofand in addition includes one or more of the following features:biologically faithful dendritic arborizations (e.g., based on digitizeddata from actual neurons), biologically accurate ion channels and theirspatially compartmentalized distributions along the dendro-somaticaxis,ligand-dependent receptors and their spatial density distriutions,dissociation constants (K_(i)) for various ligands that can bind to thereceptors (e.g., neurotransmitters and psychoactive drugs),receptor-mediated intracellular signaling through second messengersystems, and spatially compartmentalized quantitative transduction ofreceptor-mediated events to ion channel conductances. Specific valuesfor these parameters can be obtained from the published literature.

In one implemtation, the various parameter values (e.g., K_(i)) aredirectly imported into the model from online databases. Using such datasources avoids the need for extensive literature searches to collect therelevant data, and allows a more automated updated of the model as newinformation becomes available.

A schematic of the overall process for generating a biologicallyrealistic neuron is shown in FIG. 5. Although we have used the GENESISmodeling language, and although we detail the process in terms of thiscode in an exemplary embodiment, it should be understood that anyobject-oriented neural modeling language can be used in the disclosedthe methods.

Dendritic Morphology

At one time it was assumed that dendrites are simply passive conduits toconvey signals to the cell body. More recently, however, it has beenseen that dendrites are active structures that can carry out manynon-linear tasks (see, e.g., Kich & Segev (2000), Nat. Neurosci3:1171-1177). For this reason, it is important that a neural modelcontain a realistic representation of dendritic tress.

In some embodiments, neuronal morphology data can be downloaded from anonline database. For example, The Hippocampal Neuronal MorphologyArchive (found on teh world wide web at www.compneuro.org/CDROM/nmorph)is an online repository of files containing the dendritic structures ofactual cells from various brain areas (hippocampus CA1, CA3, etc). Eachfile contains a complete description of an actual cell's dendritic tree,specified using x, y, and z coordinates, diameters of dendriticsegments, branch points, and codes identifying location of thosesegments in a cell (apical vs. basal, etc). In one embodiment, thesedata are downloaded in the .swc format, a publicly available standarddeveloped and documented by MicroBrightField Inc. (Williston, Vt. USA).

Ion Channel Distribution

There is a growing appreciation that the distribution of various speciesof ion channels along the dendro-somatic axis is a key determinant ofthe manner in which neurons process signals. Along with this, there is asignificant amount of research producing detailed quantitative data inthis area (see, e.g., Migliore and Shepherd (2002), Nat. Rev.Neuroschi., 3:362-370). For example, FIG. 6 shows summary data for fourdifferent ion channels (Na⁺, I_(KA), I_(k), and Ca²⁺) for a number ofdifferent neuron types (hippocampal CA1, CA3, neocortical, etc). Panel(a) shows sodium chanel density as a function of distance from the somain different neuronal cell types. Panel (b) shows peak potassium currentas a function of distance from the soma in different neuronal celltypes. Panel (c) shows peak I_(h) cationic current density as a functionof distance from the soma in different neuronal cell types. Panel (d)shows calcium channel density as a function of distance from the soma indifferent neuronal cell types. As is clear from FIG. 6, the spatialdistribution of ion channels and ion channel currents variesdramatically between ion channel types and between neuronal cell types.

In some embodiments, ion channel distribution data can be obtained froman online database. For example, NeuronDB (found on the world wide webat senslab.med.yale.edu/senslab/NeuronDB/default.asp) lists ion channeland receptor densities for a wide variety of different cell types. For agiven cell type, a schematic form is given, as shown in FIG. 7. For eachparticular compartment (e.g., soma, [S], or proximal apical dendrites,[Dap]) ion channel and receptor densities, to the extent known, aregiven, with commentary and literature references. This information iscontinuously updated and provided via menu-driven web interfaces;downloadalbe data files are also available. Lavels of cell compartments(relative to the soma) are as follows; Dad (distal apical dendrite); Dam(middle apical dendrite); Dap (proximal apical dendrite); S (soma); AH(axon hillock); A (axon); T (axon terminal); Dbd (distal basaldendrite); Dbm (middle basal dendrite); and Dbp (proximal basaldendrite).

Implementation

In one embodiment, a program written in the C programming languageresults in a creation process that creates a GENESIS single cell modelfrom data files containing information on dendritic morphology and ionchannel densities. Specifically, the inputs to this creation program,“makeproto.c,” include (1) the aforementioned swc file describing adendritic arborization pattern, (2) a data file, “chan_dens,” listing,in tabular form, densities of all species of ion channels present atdiscrete points along the dendro-somatic axis, and (3) a data filetermed “recep_dens,” analogous to the “chan_dens” data file, thatspecifies distribution of neuroreceptors along the dendro-somatic axis.In cases lacking specific quantitative information, default values areused. Makeproto.c also outputs a “recep_compart file” that lists thedensity of each species of neuroreceptor in each of the computationalcompartments created by the program. Receptor densities are used incalculating medication effects, as described below.

The output of makeproto.c is a “p file,” which is in a format specificto GENESIS. GENESIS uses teh p file to create a functioning single-cellmodel. The p file includes lines that describe correspondingcomputational compartments of the overall model. It does so byspecifying compartment dimensions, connections with other compartments,and channel densities. Details associated with creating a p file areprovided in Bower and Beeman (1998).

When the makeproto.c program reads the swc file, it ignores all pointsthat are not branch points. It does so by neglecting “kinks” inunbranched dendritic segments. Each segment between branch points isdefined as a computational compartment. A schematic of a neuron thuscreated is shown in FIG. 8, above a drawing of the actual neuron fromwhich it was derived. The makeproto.c program places ion channels inthese compartments according to the information contained in the“chan_dens” file. If data is unavailable for a particular region, (i.e.,an area of the dendritic tree a particular distance from the soma) ionchannel or receptor density values for that region are interpolatedbetween values for the densities of the ion channels and receptors innearby regions for which they are known.

Intracellular Messaging

In many cases, when a receptor is activated (e.g., by a medication), itsultimate effects on ion channel functioning are mediated byintracellular second messenger pathways. Many sych pathways have beenidentified, and because they share common elements, there isconsiderable interaction between them. For this reason, it is difficultto predict the downstream consequences of an upstream ecen (e.g., anincrease in cAMP level resulting from activation of a G protein coupledreceptor). It is particularly difficult if there are multiple upstreamevents and significant interplay can occur. An example of such a networkof intracellular interactions is shown in FIG. 9, which summarizesreceptor→second messenger→ion channel effects in hippocampus CA1neurons.

Similations of intracellular enzymatic pathways are fairly welldeveloped. Available models generally include terms for the size of theintracellular pools of various reactants, and parameters embodyingforward and backward rate constants, which may be under enzymaticcontrol. However, these models have generally not included the far“upstream” effects of drugs at receptors, or “downstream” end effects ofsecond messengers on ion channels.

The database of quantatative cellular signaling (DOQCS) (see Sivakumaranet al. (2003), Bioinformatics, 19(3): 408-415) is an extensive onlinedatabase providing data representative of signaling cascades for a widevariety of cell types. The database (available on th world wide web atdoqcs.ncbs.res.in) includes data descriptive of a large number of secondmessengers; parameter values (e.g., for reactant pool size and rateconstants) drawn from the biological literature are included as part ofthe site.

In some embodiments of the disclosed methods, data representative ofintracellular signaling, parameters, such as that found in DOQCS, isused to model detailed intracellular interactions and their effects onion channel behaviors. The intracellular signaling data available inDOQCS is particularly convienient since much of it is both comparablewith GENESIS and freely downloaded. However, this data includes neitherdrug-induced receptor effects, nor dynamic alteration of ion channelbehaviors as discussed below.

Each element of the biologically realistic neuron model, e.g., a cellcompartment, ion channel, or synapse, is encoded as an object.Information from one object is made known to other objects via messagepassing. Additional information on GENESIS, the object-oriented languageused in disclosed methods, can be found in its online manual (availableon the world wide web at genesis-sim.org/GENESIS/Hyperdoc/Manual.html.)

Various embodiments of the disclosed methods include creating an objecttype, termed “convert,” that converts a second messenger pool size to achange in an ion channel conductance, FIG. 9, for example, shows thatincreases in arachidonic acid lead to decreases in K⁺channel activity.In this example, the size of the arachidonic acid pool (calculated bythe GENESIS object of type “pool”[genesis-sim.org/GENESIS/Hyperdoc/Manual-26html#ss26.48]) iscommunicated to a “convert” object. At each simulation time step, thisobject calculates the change in the K⁺ channel conductance from thechange in arachidonic acid level. Finally, a message reporting thischange is passed to the object representing the K⁺ channel (e.g., a“tabchannel” object) via message passing.

Receptor-level Medication Effects

For various embodiments of the disclosed methods, one of the primaryapplications is the translation of medication or neurohumoral effects tochanges in cell-level behaviors and their impact on multicellularnetwork activity (e.g., in the CA1 region of the hippocampus). Thus, thepharmacodynamics of drug-receptor interaction at the neuroreceptor levelare explicitly included in the disclosed models.

A wealth of data on the pharmacokinetic properties of many psychoactivedrugs is available in biomedical literature and, most conveniently, inonline data bases. In some embodiments, pharmacokinetic data is importeddirectly from a website maintained by the Psychoactive Drug ScreeningProgram (PDSD), an NIMH-funded initiative that provides screening ofpsychoactive composition (available on the world wide web atpdsp.cwru.edu/pdsp.htm). Both pharmacokinetic data (e.g., K_(i) values)and functional data (e.g., effects on second messenger levels, such ascAMP) are imported into the model, e.g., from an extensive onlinedatabase of these values (available on the world wide web atpdsp.cwru.edu/pdsp.asp).

Functional data is represented in the disclosed model using thefollowing formal, also used in the PDSP database (available on teh worldwide web at kidb.cwru.edu/nimh/functionP.php:

$\begin{matrix}{V = {\frac{V_{\max} \times X}{K_{act} + X} + {NSA}}} & (1)\end{matrix}$

where V=activity of the intracellular species in question;X=concentration of the composition of interest; V_(max)=maximumactivity; NSA=non-specific activity at baseline; and K_(act) is a fitingparameter. Graphically, results of functional assays appear as shown inFIG. 10.

In one embodiment, medication effects on the model neuron, are based onthe operation model of drug action summarized in FIG. 11. This approachis useful because it conceptually and mathematically separates affinity,represented by K_(A), and efficacy, as embodied by the term K_(A).Affinity is the tendency for a given molecule tp bind with a givenreceptor. (OF note, K_(A) is identical with the aforementioned K_(i). Wehave used K_(A) in teh following for the sake of consistency withKenaken (Pharmacologic Analysis od Drug-Receptor Interaction (1997), NewYork: Lippincott-Raven), a standard text.

Agonists.

FIG. 11 indicates the interdependence of agonist concentration, termed“IgA” in teh figure, and response (e.g. a change in level of anintracellular messenger) (XY plane), agonist concentration and formationof agonist-receptor complex, AR (XZ plane), and the relationship betweenAR and response, E (YZ plane).

In accordance with the operational model, response is assumed to be ahyperbolic function of agonist concentration. That is, the relationshiprepresented in the XY plane can be expressed as:

$\begin{matrix}{E_{o} = \frac{\lbrack A\rbrack E_{m}}{\lbrack A\rbrack + b}} & (2)\end{matrix}$where E_(a) is response, E_(m) is maximal response, [A] is agonistconcentration, and b is a fitting constant.The relationsip between agonist concentration and the quantity ofagonist-receptor complex, [A·R], (XZ plane) is given by

$\begin{matrix}{\frac{\left\lbrack {A \cdot R} \right\rbrack}{\left\lbrack R_{l} \right\rbrack} = \frac{\lbrack A\rbrack}{\lbrack A\rbrack + K_{A}}} & (3)\end{matrix}$Where R₁ is receptor concentration, and K_(A) is the equilibriumdissociation constant; this is the familiar Langmuir adsorptionisotherm.Combining Equations 2 and 3, and assuming K_(A)>b, the relationshipshown in the YZ plane can be expressed

$\begin{matrix}{{\frac{E_{a}}{E_{m}} = \frac{\left\lbrack {A \cdot R} \right\rbrack}{K_{E} + \left\lbrack {A \cdot R} \right\rbrack}},} & (4)\end{matrix}$

K_(E) can be calculated using

$\begin{matrix}{\left\lbrack A_{50} \right\rbrack = \frac{K_{A}}{1 + \left( {\left\lbrack R_{l} \right\rbrack/K_{E}} \right)}} & (5)\end{matrix}$where A₅₀ is the concentration of agonist that produces half maximalresponse, and can be calculated using Equation 2.

Combining Equation 4 and Equation 3 yields the following

$\begin{matrix}{\frac{E_{a}}{E_{m}} = {\frac{\left\lbrack R_{t} \right\rbrack\lbrack A\rbrack}{{K_{A}K_{E}} + {\left( {\left\lbrack R_{t} \right\rbrack + K_{E}} \right)\lbrack A\rbrack}}.}} & (6)\end{matrix}$This is the central equation of the operational model, as it definesproduction of cell-level response in terms of agonist concentration,affinity of agonist for receptor, and efficacy.Antagonists

Competitive antagonists funtion by binding to receptors. They do not, bythemselves, elicit a downstream response. Instead, they decrease thestimulation that agonists can deliver to the cell. Specifically, in thepresence of an antagonist with concentration [B] and dissociationconstant K_(B), the fractional receptor occupancy by an agonist A isgiven by the following equation:

$\begin{matrix}{\frac{\left\lbrack {A \cdot R} \right\rbrack}{\left\lbrack R_{t} \right\rbrack} = \frac{\lbrack A\rbrack}{\lbrack A\rbrack + {K_{A}\left( {1 + {\lbrack B\rbrack/K_{B}}} \right)}}} & (7)\end{matrix}$

Clearly, if a receptor is not being stimulated (by another medication ora neurotransmitter in the neurohumoral milieu) a competitive antagonistwill have no effect.

In one embodiment, a medication effects program is encoded in thecomputer language C. The output of this program is GENISIS code that isread by main_cell.g, which is the program that defines and executes theoverall cell model.

The three input files to the medication effets program are themedication characteristics file, the efficacy values file, and theneurohumoral milieu file (see FIG. 5).

Medication Characteristics

This medication characteristics input file specifies, in quantitativeterms, the agent to be screened or evaluated. The format of this file isas follows.

Line 1: concentration of agent Line 2: D1 D2 α1 α2 β1 . . . Line 3:K_(A) K_(A) K_(A) K_(A) K_(A) . . . Line 4: cAMP . . . . Line 5: V_(max)Line 6: K_(act) Line 7 NSA Line 8 PI (phosphoinositol) Line 9 V_(max)Line 10 K_(act) Line 11 NSA . . . . . .

Entries in bold are headings, those in normal typeface represent theactual data values.

Included in line 2 are all receptors at which a given agent acts. Line 3gives the equilibrium dissociation constant of the agent at eachreceptor; by convention, positive values are used for agonists andnegative values for antagonists. This constant is the concentration ofthe drug that binds to 50% of the available receptors; the lower thisconcentration, the greater the affinity the drug has for the receptor.Starting on line 4 are the changes in second messenger levels as afunction of agent concentration, as specified by the parameters V_(max),K_(act), and NSA, outlined in the preceding section. All available datafor a given composition can be included in this file.

Efficacy Values

Efficacy values are stored in a file (e.g., “K_es”) having the followingformat:

Line 1: D1 D2 α1 α2 β1 . . . Line 2: cAMP cAMP . . . Line 3: K_(E) K_(E)Line 4: PI PI Line 5: K_(E) K_(E) . . .

As described below, efficacy values (K_(E)s) are calculated in themed_fx.c program. K_(E)s in file k_es can be determined, e.g., using aK_(E) from a known interaction of the receptor and a given drug or ananalog, can be averaged over many known K_(E)s at the receptor, or canbe estimated based on molecular-level interactions.

Neurohumoral Milieu.

This file, termed neur_hum, contains data in the following format:

Line 1: neurotransmitter₁ concentration₁ Line 2: neurotransmitter₂concentration₂ Line 3: neurotransmitter₃ concentration₃

This is included as an input file both because neurohumoral stimulationcan significantly affect cell functioning, and because the effects ofantagonists at specific receptors (e.g., D2 blockers) are bestunderstood if tonic, physiological stimulation at these sites isincluded in the model.

Output of the Medication Effects Program

The medication effects program, e.g., med_fx.c, carries out thefollowing algorithm, for each computational compartment of the model.

(1) If, at the receptor in question, the agent acts as an agonist, then:

-   -   (i) If data on individual second messenger level effects are        given (i.e., there are entries on lines 4 and below in the        med_char input file described above), the effect on each second        messenger is calculated using Equation 2 (above), and    -   (ii) K_(E) is calculated using Equation 5 and Equation 2;    -   (iii) If no data on line 4 and below are available, the effect        on the messenger is calculated using Equation 6, with values for        K_(E) from the input file k_es;    -   (iv) if multiple receptors have effects on a given messenger,        individual effects are summed linearly.

(2) If, at the receptor in question, the agent acts as an antagonist,then receptor occupancy due to baseline (neurohumoral) levels ofneurotransmitter are adjusted downward using Equation 7 above. Effectsof this neurotransmitter are then calculated, using the methodologyabove.

(3) The effects of neurotransmitters in the neurohumoral milieu that actat receptors unapposed by antagonists are treated as in (1) above.

The output of the medication effects program file (e.g., med_fx.c) is afile containing levels of intracellular messengers that are expressed inGENESIS code. There are also a number of receptors that regulate ionchannels directly (i.e., not through intracellular signaling cascades).Direct effects of receptors on ion channels can be included in a processsimilar to that outlined above.

Network Model of Biologically Realistic Neurons

After implementing a biologically realistic neuron model, one thenmodels the transmission of a signal from that neuron to its neighbors.This aspect of the model can include, for example, the delays in theoutput of the neuron, either by modeling a neuron's axon by a series ofcompartments, or since the signal transmission rate along an axon isrelatively constant, by delaying the input of the axon's target by anamount proportional to the distance to that target. This aspect of themodel can also include changing the synaptic strength between neurons.

A computational model of a neural tissue sample can then be created byrepresenting the steric relationship between neurons withing the sample.This can include representing neuronal cell types present in thatsample, and the density of each of those species of neurons within thatsample.

Once the steric relationship between neurons has been represented, thenext step in constructing a model is to represent patterns ofconnectivity between neurons. This can include representing how manyneurons a particular neuron is connected to and what species thoseneurons are.

At the most schematic level, a computational model of neural tissue ismade up of a large number of independent processing units that sendmessages to each other. The nature of this model lends itself toimplementation in an object-oriented environment. In such anenvironmentm the different neuronal call types can be represented bydifferent instances of an object class. The properties of each neuronalcell type are then encapsulated as methods within the neural object.

The computational model of neural tissue can be molded in anyobject-oriented environment. An example of such an environment isGENISIS (“GEneral NEural Simulation System”). An alternative environmentfor modeling neural tissue is provided by NEURON. Both GENESIS andNEURON are widely available freeware collections of model building toolsthat can be executed on typical digital computers using a conventionaloperating system, such as LINUX.

To model neural tissue having a particular pathology, one alters thedefining characteristics of a model of normal neural tissue. Forexample, one can re-define the manner in which selected neruonscommunicate with neighboring neurons, either by altering connectivitypatterns or synaptic response characteristics. Or, one can adjust theion channel behavior or morphological characteristics of selectedneurons. This will redefine the manner in which incoming signals areintegrated to produce output signals, i.e., it will change the implicittransfer function. In an object-oriented environment, this can be done,for example, by changing the corresponding methods encapsulated by theobjects.

A psychoactive drug typically functions by either changing the implicittransfer function (e.g., by altering ion channel behaviors) of selectedneuronal cell types or changing the communication characteristics (e.g.,by changing synaptic properties) of those neurons. In particular, theimplicit transfer function and communication characteristics built intoa model of diseased neural tissue can be altered to simulate exposure ofthat tissue to a psychoactive drug. On way to determine the manner inwhich the model is to be altered pharmacologically, is to perform invitro tests of the effect of the drug on individual neurons.Alternatively, if one already knew the effect of a drug on a neuron, onecould use that knowledge and bypass the need to perform such in vitrotesting.

Once a model has been constructed, an initial exicitation is applied toselected neurons. The response has been constructed, an initialexcitation is applied to selected neurons. The respons of those neuronsto the excitation is communicated to neighboring neurons, which in turncommunicate their own responses to their own neighbors. This results ina wave of excitations that unfolds over time. The nature of theexcitation changes as the neural transfer functions and communicationsproperites change. As a result, when subjected to the same excitation,normal, diseases, and drug-exposed neural tissue provide differentoutcomes. A comparison of these different outcomes provides a basis forpredicting the effect of a particular psychoactive drug on a particulardisorder.

For example, one can alter a model of discussed neural tissue tosimulate exposure of that tissue to a drug. The altered model and theunaltered model are then subjected to the same initial stimulus (e.g., astimulus associated with a particular behavioral task). This results intwo outcomes, which will be referred to as the diseased outcome and thetreated outcome respectively. If these two outcomes are the same, onecan infer that the drug will be ineffective in treating the diseasedneural tissue.

If the treated outcome and the diseased outcome differ, it may beunclear whether the effect of the drug was beneficial or detrimental. Insuch a case, it is useful to apply the same initial stimulus to a modelof normal tissue to obtain a normal outcome. The normal outcome can thenbe compared with the treated outcome. If the treated outcome is similarto the normal outcome, one can infer that the drug will be beneficial intreating a disease modeled by the perturbation.

Drug Screening

A computational drug screening or discovery process 200, illustrated inFIG. 12, is implemented using the general computational methodologydescribed above. The process 200 includes selecting (step 202) a medicaldisorder to be studied, e.g., schizophrenia. A known drug composition isselected (step 204) as a test or potential drug composition for themedical disorder. A computational network model, in which rapidscreening of test drug compositions can be performed, is generated (step206). The effects of the test composition are then applied (step 208) tothe computational network model thus generated. The effects of the testcomposition on the network are then analyzed and evaluated (step 212).The new systems and methods of drug screening processes are discussed ingreater detail below.

Medical Disorders

Computer-based computational and neural network models provide a basisfor studying numerous psychiatric and neuropsychiatric disorders, suchas schizophrenia, Alzheimer's disorder, bipolar disorder, seizuredisorders, diffuse cerebral atrophy, and obsessive-compulsive disorder,at both microscopic and macroscopic levels. Various medical disordersand mental disorders, such as schizophrenia, Alzheimer's disease,bipolar disorder, and the like, can be modeled using the newcomputational systems and methods. Accordingly, these computationalsystems and methods can be used to screen potential drugs for thetreatment of any disorder for which behavior or physiological changescan be modeled in a network model.

Schizophrenia

Schizophrenia is a heritable disorder of the brain thought to resultfrom abnormalities that arise early in life and that disrupt the normaldevelopment of the brain. These abnormalities manifest themselves asstructural differences between a schizophrenic brain and a healthybrain. For example, schizophrenic brains tend to have larger laterialventricles and correspondingly smaller volumes of neural tissue thannormal brains. It is believed that the chemical nature of theschizophrenic brain, and in particular the manner in which it processescertain neurotransmitters, such as dopamine, is also different.

The effects of endogenous substances, such as dopamine, on normalfunctioning of the brain have been gathered by assessing how suchchemical compositions alter the physiologic functioning of neurons.Accordingly, altered physiological functioning induced by suchcompositions can be simulated in biologically realistic neurons of themodel.

Symtomatically, schizophrenia can be characterized by many functionaldeficits. Schizophrenic symptoms are divided into negative and positivecategories. Negative schizophrenic symptoms consist of behavioraldeficits such as blunting of emotions, language deficits, and lack ofmodivation. Laboratory data suggests that schizophrenics with negativesymptoms may have reduced brain activity in the prefrontal cortex.Positive schizophrenic sympoms include hallucinations, delusions, andbizzare behavior. Tests performed across a variety of modalities suggestthat delusion and hallucinations (e.g., hearing voices), are associateddysfunction in the temporal lobes, a part of the brain linked toarticulated language.

Biochemical factors are believed to underlie a number of schizophrenicand psychotic symptoms. For example, excessive production of dopamine isbelieved to play a role in the schizophrenic brain. Neurotransmittersinteracting with certain receptors mediate the transfer of chemicalinformation between neurons. It is thought that when anti-schizophrenicdrugs block dopamine receptors, the symptoms of schizophrenia arereduced. Five dopamine receptors, D1, D2, D3, D4, and D5, has beendiscovered. Their function is to bind to dopamine, which triggerschanges in the metabolic activity of the postsynaptic nerve cells.

Dopamine antagonists are drugs that block dopamine receptors. Examplesof dopamine antagonists include neuroleptics, such as chloropromazine,(THORAZINE® manufactured by GLAXOSMITHKLINE®). Various studies indicatethat controlling dopamine and dopamine receptors alleviates symptoms ofschizophrenia. THe effectiveness of clozapine, a high affinity D4receptor antagonist, has called into question the traditional assumptionthat antipsychotic medications exert their effects via action that theD2 receptor alone. However, a compound created to be a highly selectiveantagonist at the D4 receptor (L-745870) was seen to be ineffective.Antipsychotic agents like clozapine are probably effective not strictlybecause of antagonism at any one receptor, but because of theirconstellation agonisms and antagonisms at dopamine and otherneuroreceptors.

Network Models of Medical Disorders (Schizophrenia)

Referring to FIG. 13, a computational drug screening method 600 includesevaluation of the potential efficacy of a test composition as apsychiatric medication by incorporating that composition's cellularlevel effects on biologically realistic neuronal behavior in acomputational model of neural tissue having a large number ofbiologically realistic neurons. This is followed by examination of theemergent “macroscopic” behavior of the tissue. As described previously,given that such models are not exact representations of the human brain,the computational drug screening method 600 identifies test compositionsthat are candidate therapeutic compositions. Accurate assessment of testcompositions is further increased by implementation of biologicallyrealistic neurons as described above and the effects of the testcomposition on their simulated intrinsic properties.

The computational drug screening method 600 includes selecting (step602) a disorder to be studied (e.g., a psychiatric or neurologicalcondition in question such as schizophrenia). A computational networkmodel that manifests the characteristics of, in this case,schizophrenica is generated (step 604).

Once the computational network model has been generated, informationabout the known effects of the test composition or physiological neuronsis incorporated into the constituent biologically realisticcomputational neurons of the computational nework model (step 606).Details of this step are described below in connection wiht FIG. 15.

The simulated behavior of the psychiatric or neurological condition isthen observed to determine if the behavior of the computational networkmodel has been changed for the better (step 608). If the behaviorindicative of the disorder has decreased, the test composition isclassified as a candidate drug (step 610). For example, in acomputational netowrk for schizophrenia, if the network shows fewerinstances of “hallucinatory” pattern recognitions, the test compositionis classified as a candidate drug for treatment of schizophrenia (step610). Similarly, in a computational network model for Alzheimer'sdisease, improved memory recall, or a less precipitous decline in memoryperformance, would indicate that the test composition is potentially aneffective treatment of Alzheimer's disease. On the other hand, if nobeneficial behavior is observed in the network model, the effectivenessof the test composition cannot be conformed (step 612).

To generate a computational network model (step 604) manifesting thecharacteristics of the psychiatric or medical condition in question,existing models can be used or new ones can be created.

For example, one way to construct a computational model for Alzheimer'sdisease is to modify a model of a healthy brain by deleting neuralsynapses in a manner consistent with neuranatomic research findings onthe condition of neural tissue in Alzheimer's patients. This results ina computational network model that exhibits a pattern of memorydeterioration similar to that seen in patients afflicted with itsdisease. In other network models, compromised memory recall can beexhibited in a computational network model for this disorder.

In the case of a schizophrenia model, illustrated in FIG. 14, a networkmodel with “pruned” neuronal connections can be built; this is one ofthe neural abnormalities thought to underlie schizophrenia. Theresultant computational network model exhibits pattern recognitionbehavior suggestive of hallucinations.

In particular, to construct a neural computational network model thatcan simulate the cognitive functions of a schizophrenic brain, a networkmodel for a normal brain that can perform basic cognitive functions isfirst generated (step 702). Then, the network is “lesioned” (step 704)in a manner corresponding to the microanatomic changes seen in actualclinical pathological studies of schizophrenia. For example, the networkcan be lesioned by removing some of the neurons, or by alteringconnections between existing neurons. Additionally, abnormalities thatare not structural in nature, but nevertheless affect the functioning ofneurons can be applied (step 708). Such “functional lesions” caninclude, for example, degrading the functional characteristics of thecomputer-modeled biologically realistic neurons in ways analogous towhat is observed in neurons of humans afflicted with schizophrenia.Other functional lesions include simulating high dopamine levels by, forexample, incorporating the effect of dopamine on the manner in which itaffects ion channel behavior, and thus their signal integratingproperties.

If by adding structural and/or functional lesions (steps 704, 708), acomputational network model showing schizophrenic behaviors is generated(step 706), we would say a model of the disorder has been generated(step 710). Various computational network models can be implemented asdescribed in greater detail below with references to FIGS. 13-14.

By applying pathological lesions to a computational network model, onecan generate a computational network model that manifests thecharacteristics of schizophrenia (step 604). As described previously,schizophrenia effects a large range of mental activities, adverselyimpacting functions as diverse as attention, memory, languageproduction, and emotion. Schizophrenia has been associated with manydifferent brain regions, particularly the frontal areas and thehippocampus. Schizophrenia is best understood as a distributeddysfunction in which multiple brain areas pathologically interact witheach other.

Although no consensus on the etiology of schizophrenia has been reached,a number of factors are suscepted. These include neurodevelopmentaldysregulation, excessive dopamine, and various environmental stressors.In step 704 and 708 of FIG. 14, one or more pathological lesions areintroduced into the computational network model of schizophrenia in amanner consistent with these suspected etiologies.

In neurodevelopmental dysfunction of schizophrenia, excessivecortico-cortical pruning may be a factor. Thus, leasioning (step 704)can be simulated by overpruning of a normal computational network model.An example of this approach to creating a lesioned model is to use anattractor network to create a computational network model of semanticpriming, a subject is exposed to a priming word. The subject is thenasked to identify a target word that is semantically related to thepriming word. An attractor network is created in which a semantic classis operationalized as particular subset of neurons of the model. Thedetails of the model are given in Siekmeier and Hoffman, “EnhancedSemantic Priming in Schizophrenia: a Computer Model Based on ExcessivePruning of Local Connections in Association Cortex,” British Journal ofPsychiatry, vol. 160:345-350 (2002), the contents of which are hereinincorporated by this reference in their entirety.

At baseline, the semantic priming model showed semantic priming behaviorsimilar to that seen in control clinical populations. When this model isleasioned by selectively removing neural connections in way similar tothat seen in a schizophrenic brain, the percentage of correctlyidentified target words increases in a manner consistent with thatmeasured schizophrenic populations. If the priming behavior of such amodel, when exposed to a test composition, were to return to baseline,one could infer that the test compsition would be potentially effectivein alleviating symptoms of schizophrenia.

Another way to apply a lesion (step 708) is to introduce a functionaldeficit, e.g., one resulting from excessive dopamine levels in a networkmodel of a normal brain. An example of this follows.

The pre-frontal cortex and its dense dopamine activity are involved inworking memory functions. Task-related electrical activity in thepre-frontal cortex is modulated by dopamine, mainly via the D1 receptor,with dopaminergic midbrain neurons activated at the onset of workingmemory tasks and dopamine levels in the pre-frontal cortex increasingduring performance of such tasks. Blockage of dopaminergic input to thepre-frontal cortex or of dopaminergic D1 receptors in the pre-frontalcortex distributes delay-task performance. Dopamine, at physiologicallevels, appears to stablilize actual neural representations inpre-frontal cortex circuits during tasks involved in working memory,thereby rendering them robust against interfering stimuli and noise.

To mimic the dopamine-modulated ionic currents that could give rise toincreased stablility of neural representations, a network model of thepre-frontal cortex is constructed. The computational network model ofthe pre-frontal cortex includes multicompartment neurons equipped withHodgkin-Huxley-like channel kinetics that can reproduce in vitro wholecell and in vivo recordings from a pre-frontal cortex neurons.dopaminergic effects on intrinsic ionic and synaptic conductance areimplemented in the model based on the in vitro data in a variety ofways. Specifically, dopamine shifts the activation threshold of apersistent Na+ current toward hyperpolarized potentials and delays theinactivation of this current. This shift appears to contribute todopamine-induced increases in firing rate and reduced adaptationobserved in vitro., Dopamine-induced shifts of the activation thresholdcan be modeled by reducing K+ conductance; this change relies on thefact that dopamine decreases slowly inactivating K+ current inprefrontal cortex pyramidal cells. Dompamine, acting at the D1receptors, increases N-methyl-D-aspartate (NMDA)-like synaptic curentsin the prefrontal cortex. This effect can be molded as an increasedsynaptic NMDA conductance in the computational neurons. In addition,dopamine effects can be modeled by increasing glycine and g-aminobutyric acid (GABA) synaptic conductance based on dopamine inducedenhancement of GABAA (GABA, type A)-like synaptic currents in thepre-frontal cortex. Further, based on data that dopamine may heightenspontaneous firing activity of GABA neurons in pre-frontal cortex, thebackground firing rate (“noise” level) of GABAergic inputs can beincreased in the computational network model.

Consequently, dopamine strongly enhances high, delay-type activity. Thedopamine-induced changes in the biophysical properties of intrinsicionic and synaptic conductance increase stability of activatedrepresentations in pre-frontal cortex networks. Simultaneously, thedopamine-induced changes retain control over network behavior andrespond to task-related stimuli. The schizophrenic condition is simultedby exposing the model of supra-normal levels of dopamine. This causesrepresentations to become “overly” stable, resulting in responseperservation and behavioral sterotypies characteristic of schizophrenicpatients.

Other disorders can also be computationally modeled.

For many neurological and psychiatric conditions the dysfunctional brainarea may be known, but the manner in which neurobiological abnormalitiesgive rise to this dysfunction is not known. For example, it may beunclear whether the symptomatology of a particular illness arises fromdeficiencies in the GABAergic system, deficiencies in the glutamatergicsystem, abnormalities in serotonergic tone, or some combination ofthese. Thus, in the model, lesions of various magnitudes are applied toelements of each of these systems individually and in combination—thatis model parameters controlling the characteristics of the GABA system,characteristics of the glutamatergic system, and dopaminergic tone aresystematically varied and thus the “parameter space” is exhaustivelysearched.

For example, the density of GABAergic synaptic projections, can bealtered. This can be implemented by varying its density parameter in themodel through a range of values. Similarly, elements of theglutamatergic system, such as AMPA receptor density or NMDA receptordensity, can also be varied. Simultaneously, neurotransmitter levels canbe varied between subnormal and supranormal levels. As describedpreviously, the concentrations of modulatory neurotransmitters in theneurohumoral milieu are explicit inputs to the model. The manner inwhich the effects of these neurotransmitters are made known to the modelis identical to the way that medications are included. Thus, changingthe ambient serotonin level can be modeled simply by changing the valueof its concentration term.

In modeling neurological and psychiatric conditions, parameter valuesalong many dimensions can be varied—in the example given, one dimensionis related to the GABAergic system, two are related to the glutamatergicsystem, and one to the serotonergic system. With, for instance, 20gradation per parameter, this creates a large number of iterations, orparameter sets (20⁴=160,000). Using measures of model behavior thatindex pathological functioning, one can decide which of the parametersets produces behavior most like the illness in question; any modelbehavior—e.g. based on oscillatory activity or psychological task—thatdistinguishes normal from diseased behavior can be used. This, then isthe model to be used in the subsequent steps of the drug screeningprocess.

In most cases, generating a computational network model of schizophreniaincludes creating a normal computational network model (step 702) andadding structural or functional lesions to the network model (steps 704,708) to generate a computational network model of schizophrenia (step710). Referring back to FIG. 13, once the computational network model ofschizophrenia has been generated (step 604), a test composition, e.g., aknown compound, is selected and incorporated into the network (step 606)by applying information about its known effects on individualphysiological neurons to the model's constituent computational neurons.This can be carried out in a number of ways, as illustrated in FIG. 15.For example, one approach is to model (step 802) the effects of the testcompositions directly on ion channels or synaptic conductances, asillustrated in Example 3 for the drug Chlorpromazine.

A second way (step 803) is to simulated contacting neurons with drugsdrugs by way of representing their effects on cell membrane receptors,which are then transduced, via intracellular messaging, to changes inion channel behavior, this is described in detail on pages 18-28.

A third approach for applying the effects of a test composition is tomodel (step 804) the effect that the test composition is known to haveon the integration of dendritic inputs to produce an axonal output,i.e., a “transfer function.” This third method is particularly usefulwhen modeling Hopfield neurons. For example, the effects of a centralnervous system stimulatn such as methylphenidate (RITALIN®) can bemodeled by altering the transfer functions of the network's constituentneurons. Methylphenidate has a notably calming effects on hyperactivechildren and a “focusing” effect on those with attention deficitdisorders. Using this approach (step 804), methylphenidate's ability toenhance responsivity of cells can be modeled as an increase in the gainof those cell's transfer functions.

A fourth approach is to alter (step 806) model parameters associatedwith simulation of intracellular processes, such as receptor activationsand gene transcription. These processes may control receptor expression,and this method may be useful for modeling the effects ofantidepressants. A fifth approach simulates application of a testcomposition by modeling changes (step 808) in the neurotransmitterrelease properties of neurons in the network model. A sixth approachsimulates application of a test composition by changing (step 810) thesynaptic response properties of neurons in the network model. A seventhapproach simulates the manner in which a medication may change neuralconnectivity (step 812), either by enhancing neurotrophic drive,stimulating the “sprouting” of dendritic or axonal processes, or bychanging the processes by which these connections are eliminated.

In applying a test composition's physiological data to the relevantdisorder's computational model (steps 802 through 812), the informationof the test composition's effects on individual neurons or itsneuron-to-neuron (synaptic) behavior may be unknown, In such cases, itmay be necessary to first gather data experimentally by in vitroexposure of neurons to the test composition, and to measure the effectsof the test composition on individual neuron behavior and synapticresponse characteristics.

Referring back to FIG. 13, once the computational network model has beenjoined with data on the test composition's known effects (step 606), itbecomes possible to determine whether application of the drug hasmodified the behavior or functioning of the network model for the better(step 608). For example, an Alzheimer's Disease model could show adecrease in behaviors analogous to poor memory performance. In a seizuredisorder model, the network model may exhibit altered oscillatorybehavior. To this end, groups of model performance behaviors, analogousto groups of well-defined clinical performance or behavior measures, maybe employed. In particular, the gross electrical behavior of acomputational model, as measured by a simulated local field potential orsimulated electroencephalogram (EEG) can be used as an outcome measurethat is analogous to EEGs or MEGs recorded clinically, or to implantedelectrodes. The extent to which a test composition, implemented in acomputational model, decreases seizure-like model behaviors in themodel, that test composition is potentially effective at treatingepilepsy or other seizure disorders.

Also, in a schizophrenic model, the network model may exhibit alteredoscillatory behavior. For example, it has been shown using both EEG andMEG techniques, that schizophrenics have an inability to support gammafrequency (approximately 40 Hz) oscillations. Clinical experiments inwhich schizophrenic patients and controls were exposed to 20, 30, and 40Hz auditory click trains quantified this deficit, and showed that it wasunique to the 40 Hz frequency. Thus, the extent to which a testcomposition, implemented in a computational model, decreases thisoscillatory abnormality, the test composition is potentially effectiveat treating schizophrenia.

Referring to FIG. 16, a second method for using computational networkmodels to evaluate or screen potentially effective test compositionsincludes comparing the effects of test compositions with the effects ofmedications already known to be effective.

In this method (step 900), after a disorder is selected and identified(step 902), a computational model of an area of the brain known to bedysfunctional in the disorder in question, i.e., schizophrenia, isgenerated (steps 904 and 906). Two models A and B can simultaneouslygenerated. As described previously for schizophrenia, the prefrontalcortex or temporal lobes could be modeled, whereas for Alzheimer'sdisease, the basal forebrain structures could be modeled. In thisembodiment, the computational network models are not developed based onanalogies with clinically observed behaviors (as in the systems andmethods described with respect to FIGS. 14-16). Instead, theneuropathology of the brain structures involved in brains suffering fromthe disorder is modeled based on what is already known at theneuroanatomic and neurophysiologic level.

Once the computatonal network model of an area of brain disorder isgenerated (step 904), a known drug A is applied (step 908) to thecomputational network model A. As described previously, the virtualapplication of a drug compound to the constituent computational neuronsof the model can draw on: the manner in which receptor-mediated eventsneurophysiologic properties of ion channels synaptic level; data onchanges in the neurophysiologic properties of ion channels triggered byreceptor activation; or on the effects of certain classes of chemical onneuronal firing rates. Thus, the computational network model informationabout compositions that are known to be effective in treatingschizophrenia is applied (step 908) as discussed in connection with step606 of FIG. 13. For example, when the model is one of schizophrenia, theeffects of various neuroleptics (e.g., THORAZINE®, HALDOL®, ZYPREXA®)can be simulated in step 908.

Similarly, a test composition B is selected (step 910) forcomputer-based screening. The physiologic information about the testcomposition is applied (step 912) to the network model B. The resultantnetwork behavior is then examined and analyzed (step 914). After theeffects of the known effective medication or drug are applied (step908), the functioning of both network models A and B are examined bycomparing the network behaviors (step 916) of computational networkmodels A and B. Specifically, network behavior common to both the knowndrug and the test composition are noted. This can include, for example,any characteristic spatial patters of neuron activation, or anydistinctive temporal patterns in the manner in which the neural pattersof activation transition from one state to another.

To the extent that the behavior of the computational network modelsresemble each other, the test composition B can be classified as beingpotentially effective in treatment of the disorder (step 918).Conversely, to the extent that computational network models differ inbehavior, the test composition B can be classified as being ineffectivein the treatment of the disorder (step 920).

In this second drug screening method (step 900), the screening of thetest compositions is based on a direct comparison with medications knownto be effective. This second method 900 achieves this by first creatinga computational network model of dysfunctional neural tissue andsubsequently applying, to the computational network model thus created,information about drugs known to be effective in treating this disorder.The behavior of the system under the influence of the unknown testcomposition, when composed to the resultant model of the known drug,provides an indication of the test compositions potential effectivenessfor treatment.

Computer Implementations

In this embodiment described above, the computational steps of the newsystems and methods are implemented on a computer system or on one ormore networked computer systems to provide a powerful and convenientfacility for forming and testing computational network models ofbiological systems.

FIG. 17 illustrates an exemplary computer system 1000 suitable forimplementation of the new system and methods.

The computer system 1000 is illustrated as a single hardware platformincluding internal components and being linked to external components.The internal components of the computer system 1000 include a processorelement 1002 interconnected with a main memory 1004. For example,computer system 1000 can include an Intel® Pentium based processor.

The external components include a mass storage 1006. The mass storage1006 can be one or more hard disks packaged together with the processor1002 and the memory 1004. Other external components include a userinterface device 1008, which can be a monitor and keyboard, togetherwith a pointing device 1010, which can be a “mouse,” or other graphicinput devices (not illustrated). Typically, the computer system 1000 isalso linked to other local computer systems via a bus 1020, remotecomputer systems, or wide area communication networks, such as theInternet. The network link allows the computer system 1000 to share dataand processing task with other computer systems.

Several software components, which will be described in greater detailbelow, are loaded into memory during operation of this system. Thesoftware components collectively cause the computer system 1000 tofunction according to the new methods of this invention. These softwarecomponents are typically stored on the mass storage 1006. Alternatively,the software components may be stored on removable media such as floppydisks or CD-ROM media (not illustrated). A software component 1012represents an operating system (OS) responsible for managing computersystem 1000 and its network interconnections. The OS can be, e.g., ofthe Microsoft® Windows, i.e., Windows® 95, Windows® 98, Windows® NT, ora Unix operating system, such as Sun Solaris®. A software component 1014represents common languages and functions present on the computer system1000 to assist programs implementing the methods specific to thisinvention. Various programming languages that can be used to program theanalytic methods of this invention include C, C++, and the like.

The new systems and methods are programmed in mathematical softwarepackages, which allow symbolic entry of equations and high-levelspecification of processing, including special algorithms toprocedurally program individual equations or algorithms. Thecomputational models described previously can be implemented, usingfreeware packages such as GENESIS (General Neural Simulation System).GENESIS is a general purpose simulation platform developed to supportthe simulation of neural systems ranging from complex models of singleneurons to simulations of large networks made up of more abstractneuronal components. Most GENESIS applications involve realisticsimulations of biological neural systems. Other simulation software suchas NEURON, described in Hines, et al., “The NEURON simulationenvironment, Neural Comput., 9, 1179-1209 (1997), may be used. Asoftware component 1016 represents the analytical methods as programmedin a procedural language or symbolic package such as GENESIS.

Referring to FIG. 18, the software implementation of the computer system1000 may include a number of separate software components interactingwith each other. An analytic software component 1102 represents adatabase containing data for the operation of the computer system 1000.Such data generally includes, but is not necessarily limited to, resultsof prior computations, network model data, and/or clinical data. Ananalytic software component 1104 represents a data reduction andcomputational component that include one or more programs which executethe analytic methods, including the methods for testing a network model,as described in FIG. 13 and 16. Analytic software component 1106includes a user interface (UI), which provides a user of the computersystem 1000 with control and input of test network models, and,optionally, known data related to the drug screening processes. The userinterface may include a drag-and-drop interface for specifyinghypotheses to the system as shown in FIG. 1. The user interface may alsoinclude loading of network models or clinical data from the mass storagecomponents (e.g., the hard drive), from removable media (e.g., floppydisks or CD-ROM), or from a different computer system communicating withthe computer system 1000 over a network (e.g., LAN, WAN, wireless). Ananalytic software component 1108 represents the control software, alsoreferred to as UI server for controlling other software components ofthe computer system 1000.

Application of Computer Model to Drug Synthesis

As described above, a biologically realistic model of neural tissue canbe applied to modeling the effect of a drug with known characteristicson diseased neural tissue to see if a desirable outcome occurs. However,it is also possible to run this process “in reverse,” i.e., to beginwith a desirable outcome and ask what cellular level parameters must bechanged to yield that desired outcome. The answer to this questionprovides guidance as to the characteristics that are desirable for drugsthat will prove useful in the treatment of the modeled disease. Thesedesirable characteristics, in turn, inform the synthesis of novelcandidate drugs to be used for treating the disease. The method, asdescribed below, includes a process to infer a candidate drug that has aset of receptor agonisms and antagonisms efficacious for treating themodeled illness.

For a given model of a neural disorder (e.g., schizophrenia), asdescribed herein, a series of simulated drugs is generated that have arange of affinities for each neuroreceptor or ion channel between itsminimum and maximum value. In other words, the parameter space of allpossible receptor affinity values is searched vis-a-vis the varioussimulated drugs. For each simulated drug, disease outcome measures aredetermined. For example, in the model of schizophrenia, epoch frequency,semantic priming performance, are determined. The trial or trialscreating the best performance on these measures is noted. Apharmacologic agent featuring this profile of receptor activities wouldbe predicted to be most useful in treating the disease.

Stated another way, in the biologically realistic model of the diseasedtissue, there may be m parameters that can be adjusted, with each of them parameters being able to take one of n biologically realistic values.These parameters might be, for example, ion channel conductances, orsynaptic conductances. The process then identifies those values of the mparameters that result in, for example, a desired epoch frequency.Although the number of combinations of parameter values can becomelarge, it is nevertheless a finite number that can be efficientlysearched using a variety of known algorithms for searching discretevalued solution spaces of finite extent.

An analogous process can be used for single cell models. First, ahealthy (non-diseased) model neuron is evaluated on a battery ofneurophysiologic tests. This is carried out by, for example, applyingvarious inputs (e.g., a simulated fiber volley; spike train stimuli ofvarious frequencies) and quantitatively determining the neuron'sresponse properties (e.g. pattern of axonal outputs). Then, thesingle-cell model is lesioned in a manner (e.g., at the ion channel orsynaptic level) consistent with the known or suspected pathology of theneuropsychiatric illness in question. Then, the parameter space of allpossible receptor affinity values is searched, as above. The trial ortrials creating the behavior most like the normal neuron on thesemeasures is noted. A pharmacololgic agent featuring this profile ofreceptor activities would be predicted to be most useful in treating thedisease.

The following specific examples are to be construed as merelyillustrative, and not limitative of the remainder of the disclosure inany way whatsoever. Without further elaboration, it is believed that oneskilled in the art can, based on the description herein, utilize thepresent invention to its fullest extent. All publications cited hereinare hereby incorporated by reference in their entirety.

EXAMPLES Example 1 Construction of a Biologically Realistic CA1Pyramidal Cell

Using the algorithm described herein under the section entitled“Biologically Realistic Computer Model of a Neuron,”, we haveconstructed a single neuron model with spatially compartmentalizedconductances and anatomically realistic dendritic arborizations. Wecarried this out by writing a program in the C computer language thatcreates a GENESIS model based on input files containing rawneuroanatomic and neurophysiologic data, as detailed below. Input files

Dendritic arborizations. FIG. 26 shows the first 25 lines (and header)of a file in the “.swc” format, which specifies dendritic morphology ofan actual neuron. This was constructed by visualizing a rat CA1pyramidal neuron and recording the x, y, and z coordinates of pointsalong the dendritic tree. A single line in FIG. 26 contains seven piecesof information, in this order: 1.) line (segment) number; 2.) type ofsegment (1=soma, 2=axon, 3=dendrite, 4=apical dendrite, 5=fork point,6=end point); 3, to 5.) x, y, and z coordinates, respectively, of thebeginning point of segment; 6.) radius of segment; 7.)identity ofproximal segment (−1 indicates that it is the first, or most proximal,segment).

Ion channel densities. FIG. 27 is a printout of the file “chan_dens,”which lists in tabular form the densities of various species of ionchannels as a function of distance from the soma. These data weregleaned from the literature, as noted in the file (a pound sign in theleftmost position of a line indicates that it is a comment, and will notbe read by the program). Not all density values for all points along thedendritic tree are known with certainty. When it is felt safe toimterpolate between known values, an “i” is entered. If no datawhatsoever are available, a “−2” is entered; in this case, a defaultvalue is used. If not data whatsoever are available, a “−1” is entered;in this case, a default value is used. Data for a number of differentneuron types are given. The program described here uses the top eight(uncommented) lines.

Functioning neuron model. The program “makeproto.c” takes the twoaforementioned files as inputs, and creates GENESIS-readable code. Aschematic of this neuron is shown in FIG. 8. Stretches of dendrite inthe actual (biological) neuron without branch points are represented inthe computational neuron as straight segments. Ignoring “kinks” in thismanner greatly reduces the number of compartments (and computerprocessing demands) of the resultant neuron. Also, such kinks probablydo not contibute to the neuron's computational abilities in asignificant way. Thus, the number of compartments in the computationalmodel is significantly lower than the number of segments of thebiological neuron. Ion channels are distributed across the dendriticsegments based on the data of chan_dens. Each line in the GENESISreadable code output contains information fully specifying onecomputational compartment, including dimensions, ion channels present,and their densities. The resultant neuron produces a voltage trace,spiking behavior, etc. that can be compared to that of actual neurons.

Example 2 Computer Model of Normal Hippocampus

The hippocampus is divided into four subfields, CA1-4 (“CA” stands for“cornu ammonis,” another name for the hippocampus suggestive of itsresemblance to the ram's horn on the head of the Egyptian deity Ammon).

A computational model that represents tissue from the CA1 subfield ofthe hippocampus was constructed as follows.

A computer nodel representing a tissue sample from a rat brain wascreated. The virtual sample (hereafter, “the sample”) extended 154micros in the septo-temporal and transverse directions, and 634 micronsin the direction orthogonal to those two directions. The resultingsample, which extended from the stratum lacunosum-moleculare (“SL”)layer of the hippocampus to the alveus, included 400 pyramidal cells and52 interneurons. The interneurons included 16 paravalbumin (“PV”) cells,6 calbindin (“CB”) cells, 9 calretinin (“CR”) cells, and 4cholecystokinin (“CCY”) cells. This sample was then simulated by acomputational model that did not distinguish between subspecies of eachof these species of interneuron, except on the basis of their axonalprojection patterns. The dendritic morphology and ion channelsdistribution for all interneurons were assumed to be the same.

Because of constraints on computational resources, the number of cellsin the computer model was 452, only a fraction of the cells in thehippocampus. Even with this limited number of cells, twenty-four hourswere required to simulate two seconds of brain activity with ahigh-speed dual-processor personal computer. However, because thecomputer model featured many similar cell performing computationallyintensive tasks and passing results to other cells, it would have beenpossible to harness multiple computers operating in parallel. A versionof GENESIS, known as P-parallel architectures to cooperate with eachother.

Each pyramidal cell modelled using the 64 compartment pyramidal cellmodeled described by Traub, et al., “A Branching Dendritic Model of aRodent CA3 pyramidal neurone.” Journal of Physiology 1004, 481: 79-95Ppt, 1 ). The interneuron cells models were based on the model describedby Traub and Miles, “Pyramidal Cell-to-Inhibitory Cell SpokeTransduction Explicable by Active Dendritic Conductances in InhibitoryCell, ”Journal of Computational Neuroscience, 1995, 2(4):291-298. Bothmodels included a representation of the internal calcium ionconcentrations, realistic arborizations, and representation of the Na⁺,Ca⁺⁺, K⁺ _(DR), K⁺ _(AHP), K⁺ _(C), and K⁺ _(A) ion channels. Initialsegments of the axons were modeled as compartments, however axonsthemselves were modeled only as delays.

Population densities of interneuron subtypes in the model for eachhippocampal layer were calculated from known population densities.Suitable data obtained from 60 micron hippocampal slices was given byFreund and Buzsaki, “Interneurons of the Hippocampus,” Hippocampus,1996, 6(4):347-470. Densities of pyramidal cells in the stratumpyramidale of the CA1 subfield are available from Boss, et al., “On theNumbers of Neurons in Fields CA1 and CA3 of the Hippocampus of SpragueDawley and Wistar Rats,” Brain Research, 1987, 406(1-2):280-287. Withina particular stratum, the neuron distribution was modeled as randomlydistributed throughout that stratum.

The spatial distribution of synapses between axons and dendrites of thehippocampus was modeled by determining a distribution of axonalarborization associated with particular morphological classes ofneurons. Then, for each class of neurons, a distribution of synaptictargets was defined. The synaptic targets were characterized by thestratum containing the target, the targeted neuron species (i.e.pyramidal cell or interneuron), and the synaptic target area (i.e.,initial segment, soma, or dendrite). For each species, the spatialbouton densities of the axonal projection cloud were calculated. Thedistribution of synaptic targets was then used to apportion, to eachneuron of a particular species, a distribution of synapses consistentwith the spatial bouton densities for that neuron species.

Table 1, reproduced below, summarizes the computer model's assumptionsconcerning the axonal arborizations associated with neuron specieslocated in particular strata. For example, according to Table 1, themodel assumed that 45% of PV neurons in the stratum oriens (“SO”) layerof the hippocampus had basket axonal arborizations, 45% had chandelierpatterns, and the remaining 10% had bistratified axonal arborizations.CB cells were assumed to be bistratified in the SO, SR, and SL layersand radiatum-projecting in the SR layer. All CB cells in the SL layerand the SO layer were bistratified. CB cells in the SR layer were evenlysplit between bistratified and radiatum-projecting arborizations. CRcells were assumed to be interneuron projecting in the SO, SP, SL, andSR layers. SOM cells were assumed to have oriens-lacunosum molecularearborizations in the SO layer and not to exist in any other layer. CCKcells were assumed to be in the SO, SP, and SR layers and were assumedto have basket arborizations in each of those layers. PC cells wereassumed to lie only in the SP layer and to have the pyramid cellarborization in that layer.

TABLE 1 PERCENTAGE HIPPOCAMPAL AXONAL OF CELL TYPE LAYER ARBORIZATIONPOPULATION PV SO Bask 0.45 PV SO Chan 0.45 PV SO Bist 0.10 PV SP Bask0.45 PV SP Chan 0.45 PV SP Bist 0.10 CB SO Bist 1.0 CB SR Bist 0.5 CB SRRadi 0.5 CB SL Bist 1.0 CR SO Intr 1.0 CR SP Intr 1.0 CR SL Intr 1.0 CRSR Intr 1.0 SOM SO o-lm 1.0 CCK SO Bask 1.0 CCK SP Bask 1.0 CCK SR Bask1.0 PC SP Pyra 1.0

Abbreviations for axonal arborizations are as follows: “bask” means“basket,” “chan” means “chandelier (axo-axonic),” “bist” means“bistratified,” “intr” means “interneuron projecting,” “radi” means“radiatum-projecting,” and “o-lm” means “oriens-lacunosum moleculare,”and “pyra” means “pyramidal axon arborization.” Abbreviations for thefour layers for the hippocampus are as follows: “SR” means “stratumradiatum,” “SP” means “stratum pyramidale,” “SL” means “stratumlacunosum moleculare,” and “SO” means “stratum oriens.” Abbreviationsfor the neuron classes in Table 1 are as follows: “PV” means“parvalbumin,” “CB” means “CCK immunoreactive cell,” and “PC” means“pyramidal cell.”

FIG. 19 summarizes, in graphical form, the axonal and dendriticarborizations of various interneuron species. The dark circles indicatethe cell body location of each of the interneuron types; the dark linesemanating from dark circles show the orientation and laminardistribution of the dendritic tree. The hatched boxes shows the layersin which the axons of each interneuron species typically arborize. Thevertically striped boxes indicate that other interneurons, rather thanpyramidal cells, are the primary targets. Pyramidal cells are shown inoutline in the background to provide an idea of which membrane domains(somatic, proximal, or distal dendritic regions) are innervated by thevarious interneuron types. For example FIG. 19 indicates that PVinterneurons project to stratum pyramidale and the proximal area ofstratum oriens.

PV cells have been found to be predominantly either basket cells orchandelier cells. In the CA1 subfield, a small number of PV cells havebeen found to have a bistratified axon projection pattern. Thus, themodel assumed that 10% of all PV cells were bistratified, with theremainder evenly divided between basket and chandelier cells.

Certain interneurons stain as both CB and SOM cells. In the computermodel, such interneurons were modeled as SOM cells. CB cells are knownto innervate the dendrites of pyramidal cells. There are three knownsubspecies of CB cells, each of which is represented in Table 1. Exceptfor half of the CB cells in the SR, all CB cells were assumed to have abistratified axonal arborization.

The CR interneurons target the dencritic processes of otherinterneurons. Therefore, in the model, CR cells were assumed to have aninterneuron-projecting arborization, regardless of which layer they werefound in.

Within the CA1 hippocampal subfield, SOM cells are believed to be o-lmcells and to have soma within the SO layer. As shown in Table 1, themodel also made this assumption.

Virtually all CCK-immunoreactive cells are believed to have a basketaxonal aborization. As shown in Table 1, the model made the sameassumption.

Finally, PC cells are assumed to populate only the SP layer. All PCcells in the model were assumed to have pyramidal axon arborization.

Having modeled the population of cells in each layer and themorphological category to which they belong, we then modeled theconnectivity between neurons. This was done by providing the model withassumption about the targets corresponding to each axonal arborization.Assumptions concerning the targets of an axonal arborization werederived from the studies in which an axon of a particular interneurontype was labeled with an anterigrade tracer to allow visualization ofthe entire axon with all its ramifications.

Quantitative estimate of spatial synaptic density for each of theinterneuron classes were derived from studies in which axonalprojections of each of a number of different classes of interneuronswere anterigradely stained. On the basis of tissue sections, thesestudies provided an estimate of a total (linear) axonal distance pervolume and a spatial bouton (synapse) frequency, in boutons per lengthof axon.

In implementing the present model, several simplifying assumptions weremade: (i) for a given cell, bouton density was assumed to be a functionof the calcium-binding protein of the interneuron, (ii) spatial boutondensity was assumed to be invariant with septo-temporal or transversedistance from the parent cell, and (iii) for a given interneuron type,bouton density was assumed to be invariant across strata. The valuesused in the model, in units of synapses per cubic millimeter were asfollows: PV 2.12×10⁵, CB 1.13×10⁴, CR 6.36×10³, SOM 4.21×10⁵, CCK2.12×10⁵.

It is believed that CA1 axons course through the SO toward the alveus,giving off occasional collateral as they do so. Thus, in the newcomputer model, we assumed that connections between pyramid cells occuronly in the SO on the cells basal dendrites. Postynaptic targets, aspecentages, are not readily available for pyramidal cells. The modelthus assumed that the ratio of the number of connections between pyramidcells and the number of connections between pyramid cells andinterneurons was roughly proportional to their relative abundance in theSO layer.

Spatial bouton density in pyramidal axons was calculated, as describedabove in connection with interneurons to be 9.09×10⁴ synapses per cubicmillimeter.

Table 2 summarizes the model assumptions concerning connectivity betweenneurons. Each row in Table 2 corresponds to a type of axon arborization,and each column corresponds to a particular hippocampal layer. At theintersection of a row and column is a 2×3 connectivity matrix thatsummarizes the connectivity of a particular axon arborizations presentin that hippocampal layer.

TABLE 2 SO SP SR LM chan 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0bask 0 0 0 0.02 0.53 0.45 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 bist 0 00.93 0 0 0 0 0 0.93 0 0 0 0 0 0.07 0 0 0 0 0 0.07 0 0 0 intr 0 0 0 0 0 00 0 0 0 0 0.5 0 0 1.0 0 0 1.0 0 0 1.0 0 0 0.5 o-lm 0 0 0 0 0 0 0 0 0 00.03 0.86 0 0 0 0 0 0 0 0 0 0 0.11 0 radi 0 0 0 0 0 0 0 0 0.93 0 0 0 0 00 0 0 0 0 0 0.07 0 0 0 pyra 0 0 0.8 0 0 0 0 0 0 0 0 0 0 0 0.2 0 0 0 0 00 0 0 0

The first row of each connectivity matrix provides data on how theparticular axon arborization connects to pyramid cells, and the secondrow provides data on how that axon arborization connects tointerneurons. The three columns in the connectivity matrix correspond todifferent parts of the target cell to which the particular axonarborization can connect. From left to right, these are the initialsegment, the soma, and the dendrites.

For example, Table 2 shows that in the model, a bistratified axonarborization present in the SO layer was assumed to connect to adendrite on a pyramid cell with probability 0.93, and to connect to adendrite on an interneuron with probability 0.07. Basket axonarborizations were assumed to exist only in the SP layer and to connectonly to pyramid cells. Within the SP layer, the basket axon arborizationcan connect to all three parts of the pyramid cell, but only rarely toits initial segment (with probability 0.2).

The model assumptions shown in Table 2 were derived from known studiesin which axons of a particular interneuron type were labeled with ananteriograde tracer, allowing visualization of the entire axon with allof its ramifications. For example, it is known that o-lm cells projectonly to LM; that of the total number of synaptic connections they formthere, 86% are on dendrites of pyramidal cells, 3% are on pryamidal cellsoma, and 11% are on interneuron soma. Analogous data has been derivedfrom similar experiments on the other interneuron classes.

Data underlying Table 2 can be found in a variety of publications. Forexample, Katona, et al., “Postsynaptic targets ofsomatostatin-immunoreactive interneurons in the rat Hippocampus,”Neuroscience, 1999, 88(1):37-55 reports that o-lm cells project only tothe LM, and that of the total number of synaptic connections they formthere, 86% are on dendrites of pyramidal cells, 3% are on pyramidal cellsomata, and 11% are on interneuron somata. Similar data for basket cellsis found in Halasy, et al., “Synaptic Target Selectivity and Input ofGABAergic Basket and Bistatified Interneurons in the CA1 Area of the RatHippocampus,” Hippocampus, 1996, 6:306-329. Data for bistratified cellsis found in Gulyas, et al., “Pyramidal cell dendrites are the primarytargets of calbindin D28k-immunoreactive interneurons in thehippocampus,” Hippocampus, 1996, 6(5):525-34. Data fororiens-lacunosum-moleculare cells is also found in Katona, et al., citedabove. Data for pyramidal cells is found in Somogyi, P., et al.,“Salient features of synaptic organization in the cerebral cortex,”Brain Research Reviews, 1998, 26:113-135. Regardinginterneuron-projecting cells, Gulyas et al., “Interneurons containingcalretinin are specialized to control other interneurons in the rathippocampus,” J. Neurosci, 1996, 16(10):3397-411 states that CRinterneurons synapse with dendrites of other interneurons.

Because of constraints on the data structure used in implementing themodel, certain simplifying assumptions were made. For example, it isbelieved that CR interneurons synapse with dendrites of otherinterneurons, with the following caveats: (a) CB interneurons avoid PVinterneurons; (b) when CR interneurons contact other CR interneurons,they do so at both dendrites and at soma; and (c) CR interneurons formdendro-dendritic contacts (gap junctions) with one another. Thesesubleties were ignored in the present model, in part because of apaucity of experimental data for incorporating into the model.

In addition to modeling the population of neurons and their stericrelationships, were also modeled synaptic connections. In the new model,this synaptic conductance was assumed to obey the function

${{g_{syn}(t)} = {\frac{{Ag}_{\max}}{\tau_{1} - \tau_{2}}\left( {{\mathbb{e}}^{{- 1}/\tau_{1}} - {\mathbb{e}}^{{- 1}/\tau_{2}}} \right)}},\mspace{11mu}{{{for}\mspace{14mu}\tau_{1}} > \tau_{2}},$where τ₁ and τ₂ are time constants associated with the synapse and A isa normalization constant chosen such that g_(sym) reaches a maximum ofg_(max).

At the GENESIS level, neurotransmitter release into the synaptic cleftis not explicitly modeled. Rather, the synchan object is placedpostsynaptically. When a “SPIKE” message is received by the synchanobject, it activates according to parameters set in its fields.

All interneurons form GABAergic level, neutransmitter release into thesynaptic cleft is not explicitly modeled. Rather, the synchan object isplaced postsynaptically. When a “SPIKE” message is received by thesynchan object, it activates according to parameters set in its fields.

All interneurons form GABAergic synapses on their target cells.Pyramidal cells form excitatory synapses on their target cells, eitherAMPA or NMDA, with a 0.5 probability for each. Channel characteristics,in terms of values given in the foregoing equation, are presented inTable 3.

TABLE 3 E(volts) τ₁ (sec) τ₂ (sec) gmax (nS) AMPAp 0 0.0536 0.00132 8NMDAp 0 0.1440 0.00100 6 GABAp −0.075 0.0094 0.00100 3 AMPAi 0 0.05360.00132 6 NMDAi 0 0.0540 0.00130 2 GABAi −0.075 0.0094 0.00100 0.5

In the computer model, all synaptic stimulatation to cells (except forthe initial stimulation) arose from cells within the model. Because themodel represented only a very small piece of tissue, the innervation ofa typical cell in the model was considerably lower that that ofhippocampal cells in vivo. For example, actual CA pyramidal cellsreceive about 30,000 excitatory and 1,700 inhibitory inputs, which isorders of magnitude greater than the number of synaptic inputs receivedby pryamidal cells in the present model. To compensate, g_(max) for allsynaptic conductances was weighted by a dimensionless scale factor. theweight factor was obtained experimentally by gradually increasing it andre-running the computer model until pyramidal interneuron outputs showedrealistic spike waveforms. A suitable weight factor for the foregoingcomputer model was 300.

Initial Stimulus

A computer model of neural tissue, no matter how realistic, will donothing unless it receives an initial stimulus. Once subjected to astimulus, the computer model evolves through a series of states. Thisevolution depends on the various parameters defined above.

In the present example, the initial excitation was selected to bebiologically realistic, but not necessarily to contain any information.The initial excitation chosen was to set all potentials to zero and thento excite each neuron for one second using a pulse sequence having amean frequency of 15 Hz. The excitation frequency and phase were bothrandom to avoid the possibly unrealistic effect of synchronicity betweenneurons. The simulated trans-membrane potentials in each neuron in themodel were periodically recorded.

Output from Normal Model

The foregoing pattern of initial excitation produced an initialtransient period that lasted approximately 0.5 seconds. Following thistransient period, the simulated neuron potentials settled into a normalpattern of spike waveforms. Data from this transient period was excludedfrom out analysis because such data would not have been biologicallyrealistic.

The raw output provided by the model after simulation included 452time-varying voltages. Because of difficulty inherent in interpretingsuch copious data, it was useful to implant “virtual electrodes” at 16locations throughout the simulated tissue block to record local fieldpotentials averaged over a volume of tissue. The virtual electrodes werecreated by summing local field potentials at points surrounding thedesired location of the virtual electrode. The output of such a virtualelectrode was given by a sum of local field potentials at discretepoints weighted by the distances between the discrete points and thelocation of the virtual electrode:

$\Phi = {\frac{1}{4\pi\; s}{\sum\limits_{i = 1}^{n}\frac{I_{mi}}{R_{t}}}}$where Φis the field potential in volts, s is conductivity of the mediumsurrounding the neurons (in mhos), I_(mi) is the transmembrane current(in amperes) across the i^(th) neural compartment, and R_(i) is thedistance from the i^(th) neural compartment to the virtual electrode.The sum was evaluated over every computational compartment of everyneuron in the simulation.

Table 4 compares the average spike rates (in Hertz) in the model withapproximate spike rates (also in Hertz) as measured in rats engaged in avariety of tasks. These tasks include performing trials on radial 8-armmaze task (column “a”); paradoxical sleeping (column “b”), slow-wavesleeping (column “c”), and being in a “maximum” behavioral state (column“d”).

TABLE 4 model (a) (b) (c) (d) pyramidals 2.45 2 0.1-0 0.5-5 <20interneurons 25.6 15   30-100   10-40 30-100

The average spike rate for all pyramidal cells in the model was on theorder of 2.45 Hz, and for all interneurons to be 25.6 Hz. Of the 452cells in the model, 122 spiked at least once during the simulation. Thisis consistent with published results for corresponding in vivo spikerates.

The distribution of spike rates in vivo also matches the correspondingdistribution of spike rates in the model. For example, FIG. 20A shows ahistogram distribution of spike rate for pyramidal cells in the model.FIG. 20B shows a histogram of corresponding spike rates as measured invivo in a rat hippocampus. FIG. 20C and 20D are corresponding histogramsfor interneuron spike rates in both the computer model (FIG. 20C) and inthe rat (FIG. 20D).

The distribution of intervals between spikes in vivo also matches thecorresponding distribution of intervals in the model. For example, FIGS.21A-21D show the distributions of the intervals between spikes forpyramid cells (FIGS. 21A-21B) and interneurons (FIGS. 21C-21D). FIGS 21Band 21D show measured interval distributions as measured in vivo from arat hippocampus, and FIGS. 21A and 21C show modeled intervaldistributions. Comparison of these figures indicates that the intervaldistribution of the model matches that obtained in vivo.

By observing the evolution of the spatial distribution of voltagepotentials over all the neurons, it was possible to identify activity“epochs” of relative stability, separated by short transition periods.In FIG. 22B, a firing pattern for a subset of twenty neurons reveals anumber of these stable epochs. The lengths of these epochs did notappear to change for small changes in the parameters of the model.

Inspection of the spike train of a single neuron in many cases revealedepisodes of fast firing followed by episodes of markedly slower firing.Using this data, we defined a transition point as occurring when thesequential inter-spike interval ratio (ISI(t)/ISI(t+1)) was either verylarge, indicating the beginning of a fast-spiking episode, or verysmall, indicating end of such an episode. When a sufficient number(referred to as a “threshold number”) of neurons experienced atransition point within an adequately narrow time window, a transitionbetween activity epochs was said to occur. When using parameter valuesfor sequential ISI ratio, threshold, and time window of 14, 10, and 33msec, respectively, the system exhibited seven transitions, whichdefined eight activity epochs. This behavior was robust to small changesin parameter values. It was seen, to a greater or lesser degree, in alllocal field potentials measured by all the virtual electrodes in themodel. A representative trace, from a virtual electrode in the stratumpyramidale is shown in FIG. 22A.

Inspection of the data indicated that many transitions between epochswere preceded by the progressively slower firing of a number of cells.This is consistent with spike frequency adaptation (SFA), a process thatis mediated by slow Ca⁺-activated K⁺ channels. Removing the I_(AHP)channels from the constituent meurons of the hippocampus model disablesSFA. When this was done, and the model subjected to the same initialexcitation, only three activity epochs were present in the same timeperiod. A portion of the resulting firing pattern is shown in FIG. 23.

The stable states, or epochs, behave consistently with the attractorstates described in the neural network literature. As discussed earlier,in a neural network, attractors are stable fixed points of the systemthat correspond to minima of the network's energy function over allvalues of its state-space. It has been theorized that these states,occurring in biological neural tissue, represent memories. Hence, thematrix of neural connections is the storage medium for these memories.

It is therefore useful to consider whether the epochs are in factattractors. To do so, we defined a weight matrix, W, and an input vector{right arrow over (V)}_(ni).

The weight matrix was a 452×452 element array W of elements W_(ij). EachW_(ij) represented the strength of synaptic stimulation (excitatory orinhibitory) that neuron i exerted on neuron j. For example, if neuron iwere a pyramidal cell that sent x synaptic contacts to neuron j, thenW_(ij)would have been x. If neuron i were instead an interneuron, thenW_(ij) would have been −x.

The input vector {right arrow over (V)}_(ni) had the same number ofelements as there were neurons in the model (in this case 452 elements).This vector represented the normalized activation pattern of all cellsin the system during a particular epoch, n. {right arrow over (V)}_(ni)was obtained by calculating the spike rate of cell i during epoch n,then normalizing the spike rates of all cells so that all elements of{right arrow over (V)}_(ni) lay in the interval [0,1]. This wasconsistent with the assumption of rate-coded information.

The simulation was implemented using the GENESIS neural modelinglanguage and run under LINUX on a dual-processor PC as described inBower, J. M. and D. Beeman, The Book of GENESIS: Exploring RealisticNeural Models with the GEneral NEural SImulation System, 1995, SantaClara, Calif. The individual pyramidal cell and interneuron models wereported to GENESIS by Sampat and Huerta and Menschik and Finkel,respectively. Both models are available on the internet atgenesis-sim.org/BABEL/babeldirs/cells. C programs to specify cellplacement and connectivity are designed to read parameter files in theform presented in Tables 1 and 2. This allows the model to be updated asadditional neuroanatomic information becomes available. In addition,this allows the model to be selectively altered to simulate the effectsof lesioning.

Example 3 Computer Model of Schizophrenic Hippocampus

Although the neuropathological lesion of schizophrenia is not know withcertainty, postmortem studies on the hippocampi of human schizophrenicsindicate abnormalities in writing distribution of cell types,particularly interneurons, and in particular, PV interneurons. On thebasis of this published data, schizophrenia was simulated by reducingthe population of interneurons by 56%. The interneurons removed were PVinterneurons because it is believed that the hippocamus of aschizophrenic has fewer such interneurons. The resulting model will bereferred to as the “schizophrenic model.”

The remaining interneurons may show subtle changes in connectivity.These changes can be incorporated in the model by altering the data inTable 2. However, any such changes in connectivity remain unaccountedfor in the present model of the schizophrenic hippocampus, in partbecause of a paucity of suitable experimental data.

When the schizophrenic model was re-executed with the same initialstimulus as was applied to the normal model, the number of epochsobserved during the run time was reduced from eight to two. Thissuggested that the frequency of such stable epochs was associated withthe presence of schizophrenia.

Representative data following application of the stimulus to theschizophrenic model is shown in FIG. 24. Unlike the 5.5 seconds offiring data shown in FIGS. 22B and FIG. 23, FIG. 24 shows only 4.5seconds of firing data. In addition, the neurons shown in FIG. 24 areprimarily pyramidal cells, whereas the neurons shown in FIGS. 22B and 23are primarily interneurons. Nevertheless, it is clear that the firingpattern in FIG. 24 shows many fewer epochs than does the firing patternin FIG. 22B. This result indicates that the frequency of transitionsbetween epochs, or equivalently the frequency of epochs may correlatewith the presence of schizophrenia. In particular, a lower epochfrequency appears to be associated with schizophrenia. As such, a testcomposition that increase the epoch frequency potentially an effectivetreatment for schizophrenia.

Example 4 Computer Model of Chlorpromazine-Treated SchizophrenicHippocampus

If, for a given compound, receptor affinities, intracellular effectormechanisms, and ultimate effects on ion channel functioning are known,then the methods outlined above can be used. However, even if theseparameters are unknown, drug effects can be supplied, as shown in thisExample and Example 5.

The effect of chlorpromazine on neural tissue was modeled by chaning theconductances of the ion channels of each neuron. A computational modelof the manner in which the conductance of sodium and potassium ionchannels change with respect to time and voltage is given by equations4.8-4.9 in the Book of Genesis (supra) as follows:G_(Na)={right arrow over (g)}_(Na)m³h andG_(k)={right arrow over (g)}_(k)n⁴G_(Ca)={right arrow over (g)}m²h

The conductances of a sodium channel, for example, is the product of aconstant term {right arrow over (g)}_(Na), and two time-varying terms,m³ and h. A computational model of the time-variation of these terms andtheir voltage dependence is given by equations 4.11-4.13 in the Book ofGenesis:

$\begin{matrix}{\frac{\mathbb{d}m}{\mathbb{d}t} = {{\alpha_{m}(V)\left( {1 - m} \right)} - {{\beta_{m}(V)}m}}} \\{\frac{\mathbb{d}h}{\mathbb{d}t} = {{{\alpha_{h}(V)}\left( {1 - h} \right)} - {{\beta_{h}(V)}h}}} \\{\frac{\mathbb{d}n}{\mathbb{d}t} = {{{\alpha_{n}(V)}\left( {1 - n} \right)} - {{\beta_{n}(V)}n}}}\end{matrix}$

The functional form of the alpha(V) and beta(V) term are given byequations 4.21 and 4.22 in the Book of Genesis:

${\alpha_{n}(V)} = \frac{n_{\infty}(V)}{\tau_{n}(V)}$${\beta_{n}(V)} = \frac{1 - {n_{\infty}(V)}}{\tau_{n}(V)}$

In GENESIS, most ion channels are implemented using the “tabchannel”object. This object allows specifications of such fields as {right arrowover (g)}, and for specifying the α and β functions.

A computer model of neural tissue exposed to chlorpromazine was obtainedby altering the constant term of the conductance equation shown above.In particular, the schizophrenic model is altered by:

-   -   halving the value of {right arrow over (g)}_(Na) for all sodium        ion channels;    -   by multiplying the {right arrow over (g)}_(Ca) value for all        calcium ion channels by 0.71;    -   by multiplying the {right arrow over (g)}_(k) value associated        with the sustained potassium current (I_(A) or K_(A)) channel by        0.6.

For interneurons, {right arrow over (g)}_(Na) was set to 0.0349 siemens,{right arrow over (g)}_(k) was set to 0.0174 siemens, and {right arrowover (g)}_(Ca) was set to 0.0047. For pyramidal cell, {right arrow over(g)}_(Na) was set to 0.1229 siemens, {right arrow over (g)}_(k) was setto 0.0615 siemens, and {right arrow over (g)}_(Ca) was set to 0.0164siemens.

As data on the time-varying effect of chlorpromazine becomes available,the computer model can be further enhanced by altering the time-varyingterm consistent with that data.

The result of re-executing the model with the foregoing changes is shownin FIG. 25. Comparison of FIGS. 24 and 25 shows that the number ofstable epochs is increased from three in FIG. 24 to five in FIG. 25.This suggests that chlorpromazine affects the neural tissue of theschizophrenic model in a way that drives it toward normalcy.

Example 5 Model of Drugs that Act on Receptors

Not all psychiatric drugs are amenable to modeling by altering ionchannel conductances as described above. For other psychiatric drugs,only the effect on neurotransmitter receptors is known. If a drug isknown to block a particular dopamine receptor with a known affinity, onecan alter the model to simulate this effect. The method described inthis example can be brought to bear if a receptor's ultimate effects arenot mediated by intracellular messaging, or the details of suchmessaging pathways are not completely understood.

In particular, affinities of many drugs for receptors are well known.For example, Goodman and Gilmans “The Pharmacological Basis ofTherapeutics,” Tenth Edition. McGraw-Hill, 2001, list such affinites.The following equation, from Cooper, J. R. et al., “The BiochemicalBasis of Neuropharmacology,” Seventh Ed., Oxford University Press, 1996,gives the fraction of receptors occupied (r) as function ofconcentration of drug ([L]) and a dissociation constant Ki.

$r = \frac{\lbrack L\rbrack}{\lbrack L\rbrack + {Ki}}$

The Ki value, which is inversely proportional to affinity, is theconcentration needed to saturate half of the receptor. In Table 20-2(entitled “Potencies of Standard and Experimental Antipsychotic Agentsat Neurotransmitter Receptors”) of Goodman and Gilman's (supra), Ki isgiven in nanomolar (nM) quantities.

To model the effect of olanzapine on the K+channel via its antagonism ofthe D1 receptor, one first obtains olanzapine's Ki for the D1 receptor,which is 31.0. A reasonable value for the free concentration ofolanzapine (corrected for binding to plasma proteins) is 27 nM. Thisvalue is from Tran, et al, “Olanzapine (Zyprexa)—A Novel Antipsychotic,”Lippincott Williams and wilkins, 2000. The fraction of occupiedreceptors, r, is therefore is 0.47.

From Kuzhikandathil, E. V. and Oxford, G. S., “Classic D1 dopaminereceptor antagonistR-(+)-7-chloro-8-hydroxy-3-methyl-1-phenyl-2,3,4,5-tetrahydro-1H-3-benzazepinehydrochloride (SCH23390) directly inhibits G protein-coupled inwardlyrectifying potassium channels, ” Molecular Pharmacology. 62(1):119-26,2002, we see that complete D1 antagonism results in a decrease in themagnitude of the K+ current by about 30% . Similar data is available formany receptor-ion channel combinations.

Statistically, approximately 47% of D1 receptors will be blocke(antagonized) at a therapeutic olanzapine level. Therefore, the completeD1 antagonist effect described above is corrected by this amount bydecreasing the constant ({right arrow over (g)}) term in the equationfor conductance. Thus, the baseline {right arrow over (g)} is weightedby 0.30×0.47 to arrive at the {right arrow over (g)} value for theolanzapine treated neuron.

OTHER EMBODIMENTS

It is to be understood that while the invention has been described inconjunction with the detailed description thereof, the foregoingdescription is intended to illustrate and not limit the scope of theinvention, which is defined by the scope of the appended claims. Otheraspects, advantages, and modifications are within the scope of thefollowing claims.

1. A computer-implemented method for screening a receptor antagonisttest composition for potential efficacy in treatment of schizophrenia,the method comprising causing a computer to execute instructions for:(a) creating a first computer model representative of a volume ofdisease-afflicted neural tissue comprising a plurality of interconnectedmodeled neurons; wherein for each modeled neuron the modelquantitatively simulates spatially compartmentalized conductance of aplurality of ion channels and spatially compartmentalized activity of aplurality of receptors; (b) modifying the first computer model tosimulate exposure of each of the receptors R in the disease-afflictedneural tissue to the test composition by (1) defining (i) an antagonisttest composition concentration B, (ii) an antagonist test compositiondissociation constant K_(B), (iii) a naturally present agonistconcentration A, (iv) a naturally present agonist dissociation constantK_(A), (v) a receptor concentration R_(t), (vi) an equilibriumdissociation constant K_(A), and (vii) an efficacy constant K_(E); (2)calculating for each receptor R in the model a fractional occupancy ofthe receptor by naturally present agonist A in the presence ofantagonist test composition B using equation (1): $\begin{matrix}{\frac{\left\lbrack {A \cdot R} \right\rbrack}{\left\lbrack R_{t} \right\rbrack} = \frac{\lbrack A\rbrack}{\lbrack A\rbrack + {K_{A}\left( {1 + {\lbrack B\rbrack/K_{B}}} \right)}}} & (1)\end{matrix}$ wherein A·R represents a concentration of anagonist-receptor complex; (3) using the value of A·R from step (2),calculating a test composition response ratio E_(a)/E_(m) for eachreceptor in the model using equation (2): $\begin{matrix}{\frac{E_{a}}{E_{m}} = \frac{\left\lbrack {A \cdot R} \right\rbrack}{K_{E} + \left\lbrack {A \cdot R} \right\rbrack}} & (2)\end{matrix}$ wherein E_(a) is a test composition response and E_(m) isa maximal test composition response; and (4) combining results for eachreceptor in the model to provide an overall test composition efficacy;(c) providing an initial excitation to the first modified computermodel; (d) determining a first outcome indicative of a response of thefirst modified computer model to the initial excitation and the testcomposition, wherein the first outcome is representative of the efficacyof the test composition at treating the disorder; and (e) providing anoutput indicative of the efficacy of the test composition at treatingthe disorder to a physical memory, a mass storage device, a removablemedia, or a user interface device.
 2. The method of claim 1, furthercomprising: creating a second computer model representative of thevolume of disease-afflicted neural tissue comprising a plurality ofinterconnected modeled neurons each comprising a dendro-somatic axis;wherein the volume is isolated from the test composition; and whereinfor each modeled neuron the model quantitatively simulates spatiallycompartmentalized conductance of a plurality of ion channels andspatially compartmentalized activity of a plurality of receptors; andwherein each modeled neuron is generated from data describing two ormore of (1) a distribution of neuroreceptors along the dendro-somaticaxis; (2) an anatomically realistic dendritic arborization pattern, and(3) a density of ion channels present along the dendro-somatic axis;providing the initial excitation to the second computer model;determining a second outcome indicative of a response of the volume ofdisease-afflicted neural tissue represented by the second computer modelto the initial excitation; and classifying the test composition as acandidate composition for treating the disorder on the basis of adifference between the first and second outcomes.
 3. The method of claim1, further comprising: creating a second computer model representativeof the volume of disease-afflicted neural tissue comprising a pluralityof interconnected modeled neurons each comprising a dendro-somatic axis;wherein the volume is isolated from the test composition; and whereinfor each modeled neuron the model quantitatively simulates spatiallycompartmentalized conductance of a plurality of ion channels andspatially compartmentalized activity of a plurality of receptors;providing the initial excitation to the second computer model;determining a second outcome indicative of a response of the volume ofdisease-afflicted neural tissue represented by the second computer modelto the initial excitation; and classifying the test composition as acandidate composition for treating the disorder on the basis of adifference between the first and second outcomes.
 4. The method of claim1, further comprising: creating a second computer model representativeof the volume of neural tissue, wherein the neural tissue is free of thedisease and comprises a plurality of interconnected neurons, wherein foreach neuron, the model quantitatively simulates spatiallycompartmentalized conductance of a plurality of ion channels andspatially compartmentalized activity of a plurality of receptors;providing the initial excitation to the second computer model; anddetermining a second outcome indicative of a response of the volume ofdisease-free neural tissue represented by the second computer model tothe initial excitation; and classifying the test composition as acandidate composition for treating the disorder on the basis of asimilarity between the first and second outcomes.
 5. The method of claim4, further comprising: determining a first epoch frequency associatedwith the first outcome; determining a test epoch frequency associatedwith the second outcome; and classifying the test composition as acandidate composition for treatment of the disorder on the basis of arelationship between the first epoch frequency and the test epochfrequency.
 6. The method of claim 1, further comprising selecting thevolume to be a volume of a hippocampus.
 7. A computer system forscreening a test composition, the system comprising: a processor; astorage medium coupled to the processor, the storage medium encodingsoftware that, when executed, causes the processor to carry out themethod of claim
 1. 8. A computer-readable medium having encoded thereonsoftware functionally interacting with a system for screening a testcomposition for potential efficacy in treatment of a disorder, thesoftware comprising instructions for causing the system to carry out themethod of claim 1; wherein the computer-readable medium comprises amemory in a computer, a mass storage component, or a media that can beinserted or removed from a computer.
 9. The method of claim 2, furthercomprising: creating a third computer model representative of the volumeof neural tissue, wherein the neural tissue is free of the disease andcomprises a plurality of interconnected neurons, wherein for eachneuron, the model quantitatively simulates spatially compartmentalizedconductance of a plurality of ion channels and spatiallycompartmentalized activity of a plurality of receptors; providing theinitial excitation to the third computer model; and determining a thirdoutcome indicative of a response of the volume of disease-free neuraltissue represented by the third computer model to the initialexcitation; and classifying the test composition as a candidatecomposition for treating the disorder on the basis of a similaritybetween the first and third outcomes.